VAPOR  BUBBLE  GROWTH 
IN  HETEROGENEOUS  BOILING 


By 

WEN-CHIN  CHEN 


A DISSERTATION  PRESENTED  TO  THE  GRADUATE  SCHOOL 
OF  THE  UNIVERSITY  OF  FLORIDA  IN  PARTIAL  FULFILLMENT 
OF  THE  REQUIREMENTS  FOR  THE  DEGREE  OF 

DOCTOR  OF  PHILOSOPHY 

UNIVERSITY  OF  FLORIDA 


1995 


TO  MY  PARENTS,  MY  WIFE  YI-YU,  AND  MY  SON  ERIC 


ACKNOWLEDGMENTS 


I would  like  to  thank  Dr.  Klausner,  chairman  of  the 
supervisory  committee,  for  all  the  support  and  encouragement 
during  this  research.  His  enthusiasm  and  great  scholarly 
attainments  have  contributed  substantially  to  this  work. 

I also  want  to  express  my  deepest  gratitude  to  Dr.  Mei, 
cochairman  of  the  supervisory  committee,  for  his  patience, 
invaluable  guidance,  and  editing  comments  during  this  study, 
without  his  great  scholarly  achievements  this  dissertation 
would  not  have  been  possible. 

I would  like  to  express  my  appreciation  to  other 
committee  members.  Dr.  Hsieh,  Dr.  Kurzweg,  and  Dr.  Sherif. 

I also  want  to  thank  Dr.  Hongjun  Li  and  fellow  graduate 
students  and  friends  for  their  valuable  help. 

Finally,  my  deepest  gratitude  goes  to  my  parents  and  to 
my  wife,  Yi-Yu,  especially  for  her  love,  patience,  and 
understanding . 


• • • 


1X1 


TABLE  OF  CONTENTS 


ACKNOWLEDGMENTS  iii 

NOMENCLATURE  vii 

ABSTRACT  X 

CHAPTERS 

1 INTRODUCTION  1 

2 VAPOR  BUBBLE  GROWTH  IN  SATURATED  POOL  BOILING 

— PART  I.  FORMULATION 8 

2 . 1 Background 8 

2.2  Formulation 10 

2.2.1  Assumption,  Governing  Equations,  and 

Boundary  and  Initial  Condition  . . 10 

2.2.2  Non-dimensionalization  and  Coordinate 

Transformation  19 

2.3  Results  and  Discussions 27 

2.3.1  Determination  of  the  Empirical 

Constant  C and  C^ 27 

2.3.2  Comparison  of  the  Predicted  and 

Measured  Growth  Rate 35 

2 . 4 Summary 41 

3 VAPOR  BUBBLE  GROWTH  IN  SATURATED  POOL  BOILING 

— PART  II.  GROWTH  RATE  AND  THERMAL  FIELDS  . . 43 

3.1  Background  43 

3.2  Effect  of  Ja,  Fo,  k and  a on  R(r)  . . . 45 

3.3  Thermal  Field  of  the  Liquid  Microlayer  . . 51 

3.4  Thermal  Field  of  the  Solid  Heater  ....  57 

3.4.1  Effect  of  Ja 62 

3.4.2  Effect  of  Fo 66 

3.4.3  Effect  of  68 

3.4.4  Effect  of  a 71 

3 . 5 Summary  74 

4 A SIMPLIFIED  MODEL  FOR  PREDICTING  VAPOR  BUBBLE 

GROWTH  RATE  IN  HETEROGENEOUS  BOILING  77 


IV 


77 

79 

79 


4 . 1 Background  

4.2  Analytical  Formulation  

4.2.1  Vapor  Bubble  and  Liquid  Microlayer 

4.2.2  Lumped  Analysis  for  the  Solid 


Heater 82 

4.2.3  Non-dimensionalization  83 

4.3  Results  and  Discussion 85 

4.3.1  Results  for  Fo>l;  Determination  of 

the  Empirical  Constant  C2  ....  85 

4.3.2  Results  for  Fo>l;  Modified  Fourier 

Number  (Fo„) 92 

4.4  Summary  95 

5 VAPOR  BUBBLE  GROWTH  IN  SUBCOOLED  HETEROGENEOUS 

BOILING 96 

5.1  Background  96 

5.2  Theoretical  Formulation  98 

5.2.1  Vapor  Bubble  and  Liquid  Microlayer  98 

5.2.2  Subcooled  Bulk  Liquid 100 

5.2.3  Solid  Heater 103 

5.2.4  Non-dimensionalization  and  Grid 

Generation 104 

5.3  The  Asymptotic  Analysis  of  109 

5.4  Results  and  Discussion 114 

5.4.1  Validation  of  the  Model  by  the 


Comparison  of  Growth  Rate  ....  115 

5.4.2  Behavior  of  the  Heat  Transfer  Rate 

to  the  Bubble  during  Growth  . . . 119 

5.4.3  Thermal  Field  in  the  Solid  Heater 


during  the  Growth 119 

5.4.4  Thermal  Field  of  the  Bulk  Liquid  . 122 

5.4.5  Asymptotic  Behavior  of  during 

Rapid  Growth 125 

5.4.6  Velocity  Field  of  the  Bulk  Liquid  127 

5.5  Summary 129 

6 CONCLUSIONS  AND  RECOMMENDATION  FOR  FUTURE 

RESEARCH 131 

6.1  Achievements 131 

6.2  Suggestions  for  Future  Work 133 

REFERENCES 135 

BIOGRAPHICAL  SKETCH  140 


V 


NOMENCLATURE 


area  of  microlayer  wedge 

surface  area  of  vapor  bubble  dome 

c 

bubble  shape  parameter 

Cl 

microlayer  wedge  angle  parameter 

C2 

to  described  the  region  of  the  heating  solid 
influenced  by  the  microlayer  evaporation 

liquid  specific  heat 

Cps 

solid  specific  heat 

f(C) 

bubble  volume  factor 

Fo 

Fourier  number 

Fo,n 

m 

modified  Fourier  number 

h 

vertical  distance  form  the  center  of  bubble  to  the 
wall 

ti 

bubble  vertical  velocity 

Ng 

liquid  latent  heat  of  vaporization 

H 

solid  heater  thickness 

& kg 

liquid  and  solid  thermal  conductivity 

Ja 

Jacob  number 

L(r) 

the  liquid  microlayer  thickness 

Pr, 

liquid  Prandtl  number 

T evap 

the  rate  of  heat  removal  from  the  heating  surface 
due  to  microlayer  evaporation 

q// 

heat  supply  from  the  bottom  of  the  heater 

q/// 

the  internal  heat  generation  per  unit  volume 

Vi 


r 


radial  coordinate 


R' 

R(t) 

Rb(t) 


Rc(tc) 

R(t) 

t 


t 


C 


& T 

t s 


T 


subO 


U 


r 


U 


z 


dimensionless  coordinate  in  the  solid  heater 
dimensionless  coordinate  in  the  bulk  liquid 
radial  direction  in  the  spherical  moving  coordinate 
bubble  radius 

radius  of  the  liquid  microlayer  underneath  the 
bubble 

bubble  radius  for  a constant  wall  superheat 

dimensionless  bubble  growth  rate 

time 

one-g  boiling  vapor  bubble  departure  time  estimate 

vapor  bubble  growth  and  collapse  time  scale 

bulk  liquid  temperature 

liquid  and  solid  temperature 

saturated  temperature 

initial  bulk  liquid  temperature 

average  temperature  of  boiling  surface 

bulk  liquid  velocity  in  the  r direction 

bulk  liquid  velocity  in  the  z direction 

the  bubble  volume 


z coordinate  in  the  direction  normal  to  the  heating 

surface 

z dimensionless  coordinate  in  the  solid  heater 

s 

dimensionless  coordinate  in  the  bulk  liquid 
grid  parameter  (=h/R) 


Greek  symbols 

a liquid-to-solid  thermal  diffusivity  ratio  (=ce^/a^) 


Vll 


“e  ‘ “s 


^"^satO 

0 

$ 


K 

a 

e 

^bl 

Pt  & P, 

T 

^ & ri 

f & e 


liquid  and  solid  thermal  diffusivity 
extrapolated  superheat-layer  thickness 
wall  superheat 

initial  superheat  at  incipience  (=Ty-Tgg^Q) 
microlayer  wedge  angle 
velocity  potential 

liquid-to-solid  thermal  conductivity  ratio  (=k^/kg) 
surface  tension 

spherical  angle  measured  from  the  axis  of  symmetry 
liquid  kinematic  viscosity 

dimensionless  wall  superheat  (=ATgg^  (t) /AT^g^Q) 
dimensionless  temperature  of  wall  superheat 
dimensionless  temperature  of  bulk  liquid 
liquid  and  vapor  density 
dimensionless  time  (=t/tj  or  t/t^) 
computational  coordinates  in  the  bulk  liquid 
computational  coordinates  in  the  solid 


• • t 

Vlll 


Abstract  of  Dissertation  Presented  to  the  Graduate  School 
of  the  University  of  Florida  in  Partial  Fulfillment  of  the 
Requirements  for  the  Degree  of  Doctor  of  Philosophy 

VAPOR  BUBBLE  GROWTH 
IN  HETEROGENEOUS  BOILING 

By 

Wen-Chin  Chen 
May,  1995 

Chairman:  James  F.  Klausner 

Cochairman:  Renwei  Mei 

Major  Department:  Mechanical  Engineering 

A numerical  analysis  is  carried  out  to  study  vapor 
bubble  growth  in  saturated  heterogeneous  boiling.  The 
bubble  growth  rate  is  determined  by  considering  the 
simultaneous  energy  transfer  between  the  vapor  bubble,  the 
liquid  microlayer  which  resides  beneath  the  bubble,  and  the 
solid  heater.  Finite  difference  solutions  for  the 
temperature  fields  in  the  microlayer  and  heater  are  obtained 
on  expanding  coordinates  as  the  bubble  grows.  The  predicted 
bubble  growth  rate  compares  very  well  with  the  reported 
experimental  data  over  a wide  range  conditions.  Dimensional 
analysis  subsequently  reveals  four  relevant  parameters.  A 
systematic  investigation  on  the  dependence  of  the  bubble 
growth  rate  and  the  solid  temperature  on  these  parameters  is 
presented.  The  numerical  study  for  the  solid  thermal  field 
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elucidates  the  detailed  energy  transfer  beneath  a rapidly 
growing  bubble. 

A simplified  model,  based  on  heat  transfer  through  a 
wedge  shaped  liquid  microlayer  and  a lumped  thermal  analysis 
for  a solid  heater,  is  developed  for  predicting  the  vapor 
bubble  growth  rate  in  heterogeneous  pool  boiling.  A first 
order  ordinary  differential  equation  is  obtained  for  the 
bubble  growth  rate.  The  present  bubble  growth  model 
compares  well  with  the  available  experimental  data  for 
intermediate  and  moderately  wetting  fluids. 

In  subcooled  heterogeneous  boiling,  detailed 
theoretical  and  computational  analyses  are  carried  out  to 
investigate  the  vapor  bubble  growth  process.  Vapor  bubble 
growth  process  is  considered  to  be  controlled  by  the 
unsteady  energy  transfer  process  among  four  regions.  The 
unsteady  energy  equation  in  the  subcooled  liquid  is  strongly 
coupled  with  the  momentum  equation,  which  is  approximated 
with  the  Laplace  equation  by  assuming  an  inviscid  flow.  The 
bubble  radius  R(t)  is  determined  by  simultaneously  solving 
the  temperature  fields  in  the  solid  and  the  subcooled 
liquid.  Close  agreement  between  the  predicted  and  measured 
bubble  radius  requires  specification  of  the  ratio  of  the 
microlayer  base  radius  to  the  bubble  radius.  The 
computational  study  for  the  solid  and  the  subcooled  bulk 
liquid  thermal  fields  elucidates  the  energy  transfer 
process. 
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CHAPTER  1 
INTRODUCTION 


Many  industrial  process,  utilize  nucleate  boiling  heat 
transfer  due  to  its  highly  efficient  heat  transport 
mechanism.  The  requirement  for  economic  heat  exchanger 
designs  and  concern  for  safety  necessitates  both  a 
quantitative  and  qualitative  understanding  of  the  boiling 
process.  Although  intensive  fundamental  research  has  been 
conducted  over  the  past  five  decades  in  pursuit  of  the  basic 
mechanisms  governing  nucleate  boiling  heat  transfer,  many 
uncertainties  remain  unresolved.  In  general,  the  nucleate 
boiling  process  consists  of  incipience,  growth,  and 
departure  of  vapor  bubbles.  Thus,  the  understanding  of 
vapor  bubble  growth  is  essential  to  the  heat  transfer 
mechanisms  associated  with  nucleate  boiling. 

Over  the  past  four  decades,  the  rate  of  vapor  bubble 
growth  in  heterogeneous  boiling  systems  has  been  the  subject 
of  numerous  experimental  and  theoretical  investigations  due 
to  its  governing  influence  on  heat  transfer.  A number  of 
theoretical  investigations  on  bubble  growth  rate  have  been 
conducted  for  uniformly  superheated  liquids  (Dergarabedian, 
1953;  Plesset  and  Zwick,  1954;  Forster  and  Zuber,  1954;  and 
Scriven,  1959)  and  nonuniformly  superheated  liquids 
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(Griffith,  1958;  Zuber  1961;  Skinner  and  Bankoff,  1964;  Han 
and  Griffith,  1965;  Van  Stralen,  1967;  and  Mikic  et  al., 

1969) , and  various  closed-form  expressions  for  the  growth 
rate  have  been  proposed.  Extensive  experimental  data  exist 
which  demonstrate  that  such  expressions  do  not  accurately 
predict  the  vapor  bubble  growth  rate  in  heterogeneous 
boiling.  In  particular,  the  inadequacy  of  such  expressions 
have  been  clearly  shown  in  Cole  and  Shulman  (1966)  and  Van 
Stralen  et  al.  (1975)  for  high  Jacob  number,  Staniszewski 
(1959)  and  Akiyama  et  al.  (1969)  for  low  Jacob  number,  and 
Keshock  and  Siegel  (1964)  for  low  gravity. 

Labunstov  (1963)  recognized  that  the  majority  of  heat 
transfer  to  a growing  vapor  bubble  occurs  in  a liquid 
microlayer  which  resides  beneath  its  base;  he  derived  a 
vapor  bubble  growth  expression  by  assuming  a constant  wall 
temperature.  Cooper  (1983)  also  developed  a microlayer 
bubble  growth  model  by  assuming  a constant  wall  temperature. 
Dzakowic  and  Frost  (1970)  calculated  the  vapor  bubble  growth 
rate  using  a constant  wall  temperature  evaporating 
microlayer  model  in  which  the  transient  microlayer  thermal 
field  is  accounted  for.  Srinivas  and  Kumar  (1984)  proposed 
a microlayer  bubble  growth  model  which  accounts  for  the 
bubble  shape  as  well  as  natural  convection  and  assumes  a 
constant  wall  temperature.  However,  these  models  do  not 
agree  well  with  available  experimental  data  (Cole  and 
Shulman,  1966;  Van  Stralen  et  al.,  1975;  Staniszewski,  1959; 
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Akiyama  et  al.,  1969;  and  Keshock  and  Siegel,  1964).  Cooper 
and  Vijuk  (1970)  and  Van  Stralen  et  al.  (1975)  have 
considered  the  heat  transfer  to  a growing  vapor  bubble  from 
both  the  microlayer  and  the  liquid  surrounding  the  vapor 
dome.  Asymptotic  expressions  for  the  vapor  bubble  growth 
rate  were  derived  by  combining  Rayleigh's  solution  for  the 
initial  stage  of  growth  with  the  heat  transfer  controlled 
solution.  The  growth  model  of  Van  Stralen  et  al.  (1975) 
appears  to  agree  better  with  the  high  Jacob  number  data  than 
that  of  Cooper  and  Vijuk  (1970).  However,  its  usefulness  is 
limited  due  to  the  required  specification  of  a dimensionless 
bubble  growth  parameter  for  which  little  guidance  was 
provided  on  how  it  should  be  calculated.  Fyodorov  and 
Klimenko  (1989)  recently  proposed  a microlayer  bubble  growth 
rate  model  in  which  spatial  variations  in  the  wall 
temperature  were  included,  but  the  temporal  variations  were 
ignored.  Their  growth  rate  model  predicts  that  the  bubble 
radius  is  proportional  to  the  square  root  of  time, 

However,  the  data  of  Akiyama  et  al.  (1969)  clearly 
demonstrate  that  the  exponent  on  t can  be  as  small  as  1/5  at 
very  low  Jacob  number. 

In  their  model  for  predicting  the  growth  and  detachment 
of  a vapor  bubble  attached  to  a heating  surface,  Lee  and 
Nydhal  (1989)  incorporated  the  effect  of  hydrodynamics 
surrounding  the  dome  of  the  growing  bubble.  They  assumed  a 
constant  wall  temperature,  modelled  the  bubble  as  consisting 
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of  a wedge-shaped  inicrolayer  and  a hemisphere,  and  carried 
out  detailed  numerical  computations  for  the  momentum  and 
energy  equations  for  the  liquid  phase  surrounding  the  vapor 
bubble.  However,  their  predictions  for  the  bubble  growth 
rate  and  the  detachment  radius  are  not  satisfactory  when 
compared  against  experimental  data. 

Despite  the  importance  of  predicting  the  bubble  growth 
rate  in  heterogeneous  boiling  and  continuous  research 
efforts  covering  the  subject,  no  satisfactory  theory  has 
been  formulated  to  cover  a wide  range  of  conditions.  In  the 
recent  works  for  predicting  the  vapor  bubble  detachment 
diameter  for  pool  boiling  (Zeng  et  al.,1993)  and  that  for 
flow  boiling  (Zeng  et  al.,  1993),  a successful  model  was 
developed  based  on  first  principles  of  mechanics.  The 
predictions  for  the  detachment  radius  agree  quite  well  with 
the  data  over  a wide  range  of  conditions;  however,  the 
agreement  was  obtained  only  after  using  the  actual  bubble 
growth  rate  reported  for  each  experiment.  Clearly,  without 
a reliable,  physically  sound  bubble  growth  rate  model,  the 
utility  of  the  bubble  detachment  models  (Zeng  et  al.,  1993a 
and  1993b)  is  limited. 

Vapor  bubble  growth  and  collapse  in  subcooled  boiling 
were  experimentally  studied  by  Gunther  and  Kreith  (1950), 
Ellion  (1954) , and  Cochran  and  Aydelott  (1966) . Theoretical 
studies  were  investigated  by  Plesset  and  Zwick  (1954), 
Bankoff  and  Mikesel  (1958),  and  Zuber  (1961).  In  their 


5 


model  for  predicting  bubble  growth  and  collapse  in  highly 
subcooled  nucleate  boiling,  Bankoff  and  Mikesell  (1958) 
considered  the  effects  of  liquid  inertia  only.  Plesset  and 
Zwick  (1954)  studied  vapor  bubble  collapse  in  a uniformly 
subcooled  liquid.  They  suggested  that  the  thermal  effects 
may  be  unimportant  for  the  collapse  rate.  From  Ellion's 
experimental  data,  Zuber  (1961)  considered  bubble  growth  and 
collapse  separately.  Zuber  extended  the  Bosnjakovic- Jakob 
theory  for  bubble  growth  and  used  Rayleigh's  solution  for 
bubble  collapse.  Rayleigh's  equation  is  only  applicable 
when  heat  transfer  is  negligible. 

In  the  present  work,  detailed  numerical  analyses  are 
carried  out  to  study  vapor  bubble  growth  in  saturated  and 
subcooled  heterogeneous  boiling.  In  the  study  of  saturated 
heterogeneous  pool  boiling,  the  energy  transfer  process  is 
considered  to  be  among  three  phases:  the  vapor  bubble,  the 
liquid  microlayer,  and  the  solid  heater.  The  unsteady 
energy  equations  in  the  solid  heater  and  the  liquid 
microlayer  are  solved  by  using  finite  difference  method. 
Finite  difference  solutions  for  the  temperature  fields  are 
obtained  on  expanding  coordinates  as  the  bubble  grows.  The 
bubble  growth  rate  is  obtained  simultaneously  with  the 
thermal  fields  at  every  time  step.  The  predicted  bubble 
radius  compares  very  well  with  the  reported  experimental 
data  over  a wide  range  of  conditions.  The  current  analysis, 
together  with  the  experimental  bubble  growth  rate  data. 
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indicates  that  the  liquid  Prandtl  number  has  little  effect 
on  the  bubble  growth  rate.  Through  dimensional  analysis, 
four  relevant  dimensionless  parameters  have  been  found. 

They  are  Jacob  number  Ja,  Fourier  number  Fo,  the  liquid-to- 
solid  conductivity  ratio  k,  and  the  liquid-to-solid  thermal 
diffusivity  ratio  a.  Systematic  investigation  on  the 
dependence  of  the  bubble  growth  rate  and  the  solid 
temperature  is  conducted  to  assist  in  the  basic 
understanding  of  the  bubble  growth  mechanism.  Although  the 
above  two-dimensional  axisymmetric  unsteady  bubble  growth 
model  is  accurate  for  predicting  bubble  growth  rate,  it  is 
impractical  to  use  as  a design  tool.  Thus,  a simplified 
model,  based  on  heat  transfer  through  a wedge  shaped  liquid 
microlayer  and  a lumped  thermal  analysis  for  a solid  heater 
is  developed  for  predicting  the  vapor  bubble  growth  rate.  A 
first  order  ordinary  differential  equation  is  obtained  for 
predicting  the  bubble  growth  rate.  The  present  bubble 
growth  model  compares  well  with  the  available  experimental 
data.  In  the  study  of  subcooled  heterogeneous  pool  boiling, 
the  energy  transfer  process  is  considered  to  be  among  four 
regions:  a vapor  bubble  whose  shape  is  assumed  to  be  a 
truncated  sphere,  a liquid  microlayer,  a subcooled  bulk 
liquid,  and  a solid  heater.  Due  to  the  strong  effect  of 
subcooling  and  rapid  bubble  growth  and  collapse,  the 
unsteady  energy  equation  in  the  subcooled  liquid  is  coupled 
with  the  momentum  equation;  the  momentum  equation  is 
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approximated  by  the  Laplace  equation.  A grid  generation 
technique  has  been  used  in  solving  the  temperature  field  of 
the  subcooled  bulk  liquid.  To  ensure  the  accurancy  of  the 
present  numerical  results  in  the  subcooled  bulk  liquid,  an 
asymptotic  solution  for  small  time  has  been  obtained.  The 
numerical  results  in  the  subcooled  liquid  compare  well  with 
the  asymptotic  solution.  The  computational  study  for  the 
subcooled  liquid  and  the  solid  heater  elucidates  the  energy 
transfer.  In  the  computation  two  parameters,  <p  and  ju^,  need 
to  be  specified  in  order  to  obtain  good  agreement  between 
the  predicted  and  measured  bubble  radius.  The  numerical 
results  show  that  only  one  combination  for  these  two 
parameters  is  possible  to  obtain  the  best  fit. 

This  work  is  divided  into  six  chapters.  The 
formulation  for  vapor  bubble  growth  in  saturated  pool 
boiling  is  given  in  Chapter  2.  In  Chapter  3,  systematic 
investigations  on  the  growth  rate  and  thermal  fields  are 
considered.  A simplified  model  for  predicting  the  vapor 
bubble  growth  rate  is  developed  in  Chapter  4 . Vapor  bubble 
growth  in  subcooled  heterogeneous  boiling  is  studied  in 
chapter  5.  Finally,  conclusions  and  recommendations  for 
future  work  are  discussed  in  Chapter  6 . 


CHAPTER  2 

VAPOR  BUBBLE  GROWTH  IN  SATURATED  POOL  BOILING 

— PART  I.  FORMULATION 

2 . 1 Background 

In  the  previous  work  considering  the  growth  rate  of  a 
bubble  attached  to  a heating  surface  (heterogeneous 
boiling),  the  simultaneous  spatial  and  temporal  changes  of 
the  heating  surface  temperature  were  not  considered. 

However,  the  temperature  distribution  in  the  neighborhood  of 
a heater  wall  has  long  been  believed  by  Westwater  (1963)  to 
be  important  to  the  mechanics  of  bubble  growth.  The  growth 
of  a vapor  bubble  requires  a certain  amount  of  energy  from 
the  heating  surface  to  vaporize  the  surrounding  liquid.  In 
saturated  boiling,  this  energy  must  come  from  the  microlayer 
which  resides  beneath  the  base  of  the  bubble  and  in  turn 
withdraws  energy  from  the  solid  heater  (Moore  and  Mesler, 
1961).  Typically,  the  growth  and  departure  of  a vapor 
bubble  involves  a very  short  period  of  time  (Hospeti  and 
Mesler,  1969).  Although  energy  is  continuously  supplied  to 
the  solid,  the  time  scale  required  for  the  solid  to  adjust 
its  temperature  into  a uniform  distribution  for  most  boiling 
systems  is  much  longer  than  the  bubble  departure  time. 
Consequently,  the  local  solid  temperature  directly  beneath 
the  microlayer  will  decrease,  both  in  space  and  in  time. 
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because  the  local  rate  of  energy  removal  due  to  the  bubble 
growth  is  typically  much  larger  than  that  added.  The  rapid 
decrease  of  the  heating  surface  temperature  directly  beneath 
a growing  vapor  bubble  has  been  experimentally  demonstrated 
by  Moore  and  Mesler(1961)  and  Hospeti  and  Mesler  (1969).  As 
the  local  solid  wall  temperature  decreases  with  time  during 
the  bubble  growth  period,  the  temperature  gradient  across 
the  microlayer  is  reduced  and  in  turn  directly  reduces  the 
rate  of  energy  transferred  to  the  bubble.  Therefore,  the 
actual  bubble  growth  rate  should  be  smaller  than  that 
predicted  by  assuming  a constant  wall  temperature.  This 
prognosis  is  in  qualitative  agreement  with  the  existing  data 
which  show  that  the  growth  rate  R(t)  is  less  than  or  equal 
to  in  which  /3  is  a constant. 

In  this  saturated  pool  boiling  investigation,  the 
bubble  growth  process  is  considered  to  be  controlled  by  the 
simultaneous  unsteady  energy  transfer  process  among  three 
phases:  the  vapor  bubble,  the  liquid  microlayer,  and  the 
solid  heater.  The  bubble  is  assumed  to  have  a uniform 
temperature  at  saturated  conditions,  Tga^of  "th©  liquid 
microlayer  is  assumed  to  be  very  thin;  and  the  temperature 
field  in  the  solid  is  assumed  to  be  axisymmetric  with 
respect  to  (w.r.t.)  the  axis  drawn  through  the  center  of  the 
bubble.  In  section  2.2.1,  the  pertinent  assumptions  are 
made  and  the  equations,  the  boundary  conditions,  and  the 
initial  conditions  governing  the  unsteady  energy  transfer 
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from  the  solid  layer  through  the  microlayer  to  the  bubble 
are  presented.  In  section  2.2.2,  non-dimensionalization  and 
coordinate  transformations  are  presented,  which  lead  to 
necessary  simplifications  of  the  governing  equations.  The 
solution  procedure  for  the  microlayer  and  solid  temperature 
fields  and  the  bubble  radius  is  also  described.  The 
parameters  characterizing  the  bubble  shape  and  the 
microlayer  wedge  angle  during  the  growth  are  determined  by 
matching  the  predicted  growth  rate  with  the  existing 
experimental  data.  In  section  2.3  The  predicted  bubble 
growth  rate  is  compared  with  the  available  experimental  data 
over  a wide  range  of  boiling  conditions  ; very  good 
agreement  is  obtained.  Notably,  it  is  found  that  the 
Prandtl  number  has  little  or  no  effect  on  the  bubble  growth 
rate  for  the  experimental  conditions  considered.  To  assist 
the  basic  understanding  of  bubble  growth  in  heterogeneous 
boiling,  the  results  of  a systematic  study  considering  the 
dependence  of  the  bubble  growth  rate  and  the  thermal  fields 
of  the  liquid  microlayer  and  the  solid  heater  on  the 
dimensionless  parameters  are  presented  in  Chapter  3 of  this 
investigation . 


2 . 2 Formulation 

2.2.1  Assumptions,  Governing  Equations,  and  Boundary  and 
Initial  Conditions 

2. 2. 1.1  On  the  vapor  bubble 


11 


A rigorous  description  of  the  vapor  bubble  growth  and 
heat  transfer  processes  among  three  phases  requires  a 
complete  account  for  the  hydrodynamics  around  the  rapidly 
growing  bubble  in  addition  to  the  complex  thermal  energy 
transfer.  The  most  difficult  part,  perhaps,  is  the  exact 
determination  of  the  bubble  shape  which  is  a moving 
interface  between  the  vapor  and  liquid.  The  numerical 
analysis  by  Lee  and  Nydahl  (1989)  still  relied  on  an  assumed 
shape  for  the  bubble,  although  the  hydrodynamics  based  on 
the  assumed  bubble  shape  is  properly  accounted  for.  The 
discussion  in  the  Introduction  suggests  that  the  thermal 
coupling  between  the  liquid  microlayer  and  the  solid  heater 
is  perhaps  more  important  than  hydrodynamics.  Therefore,  in 
order  to  simplify  the  problem,  the  hydrodynamics  of  the 
liquid  motion  induced  by  the  growing  bubble  is  not  pursued 
in  this  analysis.  The  determination  of  the  bubble  shape 
thus  relies  on  empirical  evidence.  Furthermore,  it  is 
assumed  that  energy  transfer  to  the  bubble  occurs  only 
through  the  microlayer,  and  energy  transfer  to  the  bubble 
dome  is  small.  Thus,  the  current  analysis  is  only 
applicable  to  saturated  boiling. 

In  this  study,  the  liquid  microlayer  between  the  vapor 
bubble  and  the  solid  heating  surface  is  assumed  to  have  a 
simple  wedge  shape  with  an  angle  <p<l.  The  interferometry 
measurements  of  Koffman  and  Plesset  (1983)  demonstrate  that 
a wedge  shaped  microlayer  is  a good  approximation.  There 
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exists  ample  experimental  evidence  (Van  Stralen  et  al.; 
1975;  Akiyama  et  al.,  1969;  and  Moore  and  Mesler,  1961) 
suggesting  that,  as  the  bubble  grows,  the  bubble  dome  has 
approximately  a spherical  shape  with  a radius  R(t).  Thus, 
the  bubble  is  modelled  as  a truncated  sphere  on  a wedge 
shaped  microlayer  as  depicked  in  Figure  2.1.  Using 
cylindrical  coordinates,  the  local  microlayer  thickness  is 
denoted  by  L(r).  The  radius  of  the  wedge-shaped  interface 
is  denoted  by  Rb(t),  which  is  typically  not  equal  to  R(t). 
Let 


and  the  vapor  bubble  volume  Vb(t)  can  be  expressed 
approximately  as 

nCO  = ^R\t)f(C) 

with 


(2.1) 


(2.2) 


/((7)  « 1 --(1  - Vl  - <y-f  + -(1  -\/l  - for  4><1. 

4 4 (2.3) 

In  the  limit  C-»l,  the  bubble  is  a hemisphere  and 
V^(t)  =—uR^(t).  In  the  limit  C^O,  the  bubble  approaches  a 
sphere  and  F^(  t)^— i?^(  t)  . The  dependence  of  C on  other 
parameters  will  be  determined  from  measurements  of  the 
bubble  shape  reported  under  various  conditions . In  this 
study,  C is  assumed  to  be  independent  of  time  for 
simplicity. 


Figure  2. 1 Sketch  for  the  growing  bubble,  microlayer  and  the  heating  solid. 
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An  energy  balance  for  the  growing  bubble  at  the 
microlayer  liquid-vapor  interface  depicted  in  Fig.  2.1  is 
described  as 


dt 


-k 


dZ  iz=L(r) 


dA,  for  (j)<l 


(2.4) 


where  is  the  vapor  density,  h^g  the  latent  heat,  the 
liquid  thermal  conductivity,  the  liquid  temperature,  and 
A the  surface  area  of  the  wedge  at  the  liquid-vapor 
interface.  With  axisymmetry  and  <p<l,  dA-2rrrdr  and  A~  nnl{t) 
at  a given  instant.  Equation  (2.4)  simply  states  that  the 
energy  conducted  from  the  liquid  to  the  bubble  through  the 
microlayer  is  used  to  vaporize  the  surrounding  liquid  and 
thus  expands  the  bubble.  The  initial  condition  for  R(t) 
will  be  discussed  later. 

2. 2. 1.2  On  the  microlaver 

The  microlayer  is  assumed  to  be  a wedge  centered  at  r=0 
with  local  thickness  L(r).  Because  the  hydrodynamics  are 
not  considered,  the  wedge  angle  <j>  cannot  be  determined  as 
part  of  the  solution.  In  Cooper  and  Lloyd  (1969),  the  angle 
(p  was  related  to  the  viscous  diffusion  length  of  the  liquid 
as  i?^(  t)tancp  = t in  which  is  the  kinematic  viscosity  of 
the  liquid.  Hence  for  small  (p,  one  has 

_ 

R,(t)  ■ CB(ty  (2.5) 

Cooper  and  Lloyd  (1969)  estimated  C^  to  be  within  0.3  to  1.0 
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for  the  conditions  under  their  consideration.  It  is 
postulated  here  that  the  behavior  of  (p  is  insensitive  to 
whether  or  not  the  wall  temperature  is  constant.  Thus  the 
above  form  of  equation  (2.5)  for  (p  is  adopted.  The  constant 
Cl  will  be  determined  from  matching  the  predicted  growth 
rate  R(t)  with  experimentally  measured  growth  rate.  The  two 
dimensionless  constants,  C and  C^,  constitute  the  only 
empiricism  of  the  present  analysis.  They  arise  because  the 
hydrodynamics  around  the  bubble  are  not  taken  into  account 
in  the  analysis  which  follows. 

For  both  the  microlayer  and  the  solid,  the  temperature 
fields  are  assumed  to  be  axisymmetric . The  energy  balance 
for  the  liquid  within  the  microlayer  is 

dT^  \ I d \ 

dt  rdr  dr  I dz^ 

where  is  the  liquid  thermal  diffusivity.  The  boundary 
conditions  are  given  as  follows 


T = T 

■‘sat 

at  z=L(r), 

(2.7) 

II 

II 

at  z = 0, 

( 2 • 8H&b ) 

CD 

II 

o 

at  r = 0, 

(2.9) 

o 

II 

at  r=Rlt). 

(2.10) 

In  the  above,  Tg  is  the  solid  temperature.  Equations 


for  Q<r<RXt),  Q<z<L{r) 


(2.6^ 
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(2.8a&b)  simply  ensure  the  continuity  of  the  temperature  and 
heat  flux  on  the  liquid-solid  interface.  Equation  (2.10)  is 
only  an  approximation  in  which  radial  heat  conduction  at  the 
edge  of  the  microlayer  is  neglected.  Equation  (2.9) 
guarantees  that  is  finite  as  r-»0.  The  initial  condition 
for  is  discussed  later. 

2. 2. 1.3  On  the  solid  heater 

The  energy  equation  for  the  solid  is 


dT\ 

dt 


for0<r<<xi,  -H<z<0 


(2.11) 


where  Og  and  kg  are  the  thermal  diffusivity  and  conductivity 
of  the  solid,  and  is  the  volumetric  heating  and  H is  the 
thickness  of  the  heater.  The  boundary  conditions  are 


II 

o 

at  r=0 

and 

(2.12a&b) 

II 

1 

at  z=  -H, 

(2.13) 

dT.  dT 

T = Tf  and  L ^ = k ^ 
^ ^ ^ dz  ^ dz 

at  z=0, 

Q<r<R,{t), 

(2.14a&b) 

dT 

T = Tf  and  ^ - 0 

^ ^ dz 

at  z=0 

r>Ri,{t). 

(2.15a&b) 

In  the  above,  two  possible  mechanisms  for  heat  supply  are 
incorporated:  a heat  supply  on  the  lower  surface  at  z=-H, 

, or  a volummetric  heat  supply  (on  average)  within  the 
solid  layer,  q^'^ . For  r>Rb(t)  at  z=0,  the  heat  transfer  due 
to  natural  convection  of  the  liquid  is  neglected.  Whether 
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or  not  the  heat  removal  from  the  solid  due  to  natural 
convection  is  included,  it  has  virtually  no  influence  on  the 
solution  because  it  is  so  small  compared  with  latent  heat 
removed  due  to  bubble  growth.  The  condition  at  r-»«>  reflects 
the  uniformity  of  the  temperature  far  away  from  the  bubble. 
Because  of  the  continuous  energy  supply,  it  is  not  quit 
correct  to  specify,  for  all  times,  Tg=Tgat+ATgato  r-x»  in 
which  ATgato  is  the  initial  superheat  at  incipience.  If 
there  exists  strong  thermal  interaction  among  adjacent 
nucleation  sites,  as  might  occur  in  the  coalesced  bubble 
regime,  the  present  analysis  needs  to  be  modified.  Hence, 
the  present  analysis  is  restricted,  strictly  speaking,  to 
the  isolated  bubble  regime. 

2. 2. 1.4  On  the  initial  conditions 

It  must  be  noted  that  the  above  formulation  requires  a 
well  established  microlayer  at  the  beginning  of  the 
solution.  This  implies  that  the  solution  must  start  at  some 
finite  but  small  time  in  which  t<j  is  the  depart  time  of 

the  bubble  and  is  the  time  scale  of  interest  for  the  present 
problem.  It  is  assumed  here  that  at  t=to  the  microlayer  has 
formed  and  the  bubble  has  a radius  Rq.  But  since  Rq  is 
small,  the  bubble  has  withdrawn  little  energy  from  the  solid 
at  t=to.  Thus,  the  solid  temperature  Tg  is  practically 
constant  everywhere. 


T = T + 

s sat 


AT 


satO 


at  t=  tj. 


(2.16) 
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This  assumption  is  justified  because  the  initial  stage  of 
the  growth  is  governed  by  hydrodynamics  considerations  as  is 
the  case  for  Rayleigh's  solution.  The  present  study  focuses 
on  the  coupling  of  the  solid  temperature  with  that  of  the 
liquid  and  the  rate  of  bubble  growth.  However,  this 
coupling  becomes  important  only  after  a significant  amount 
of  energy  is  transferred  from  the  solid  to  the  bubble 
through  the  microlayer.  Therefore,  a reasonable  initial 
condition  for  1i  is 


T,(t=t^,r<R^,z)  = + A 7;,,^ 


1 


Ur) 


3t  t— 


(2.17) 


Consistent  with  equations  (2.16)  and  (2.17),  the  initial 

bubble  radius  at  t=t(,  can  be  obtained  by  assuming  a constant 

wall  temperature.  This  is  appropriate  because  the  solid 

temperature  is  not  affected  yet  by  the  growth  of  the  bubble 

at  very  small  tg.  Integrating  equation  (2.4)  by  assuming  a 

dT 

time  independent  , one  obtains 


Rkt)  = 


Cy-  Ja 


m c^prf 


1/2 


above 


(2.18) 


for  the  growth  rate  with  a constant  wall  temperature.  In  the 


and 


/a  - ^pe^^satO 

pA, 


Pr  = 


CL. 


(2.19) 


(2.20) 


are,  respectively,  the  Jacob  number  and  Prandtl  number,  and 
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Cp^  is  the  specific  heat  of  the  liquid.  Hence  Ro  can  be  found 
as  Ro=Rc(to).  In  equation  (2.18),  the  relation  between  <t>  and 
R(t)  given  by  equation  (2.5)  has  been  used.  Equations 
(2.16-2.18)  complete  the  specification  of  the  initial 
conditions.  Although  the  actual  initial  growth  is  governed 
by  hydrodynamics,  this  assumed  initial  radius  has  little 
effect  on  the  later  growth  as  long  as  t(,<td.  It  is  noted 
that  the  liquid  Prandtl  number  Pr^  results  from  Cooper  and 
Lloyd's  model  for  the  microlayer  wedge  angle  0 given  by 
equation  (2.5),  and  using  equation  (2.18)  <p  may  be  expressed 
as 

_ f(C)  (tPr, 

C^Ja  (2.21) 

2.2.2  Non-dimensionalization  and  Coordinate  Transformation 
2. 2. 2.1  Non-dimensionalization  and  further  simplification 
To  gain  an  understanding  of  the  important  physics  and 
to  efficiently  solve  the  energy  equations,  the  following 
dimensionless  variables  are  introduced, 

T = tlt^,  z^^zIUf),  z=z!H, 

R(T)-R(t)IRlQ,  9,=(7;-rj/Ar„„,  <2-22) 

where  t^  is  the  bubble  departure  time  scale.  An  accurate 
determination  of  t^  is  not  necessary  because  only  the  time 
scale  is  relevant.  The  length  scale  in  the  z-direction  for 
the  solid  is  the  thickness  of  the  heater.  However,  this  is 
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The 


not  the  correct  length  scale  if  H is  very  large  or  the  solid 
conductivity  is  very  small.  Thus  another  transformation  is 
expected  for  5^  in  carrying  out  the  computation, 
utilization  of  H is  sufficient  for  the  purpose  of 
dimensional  analysis. 

Upon  substitution  of  equation  (2.22)  into  equation 
(2.4),  the  normalized  bubble  growth  rate  is  governed  by 


dR^{T')  _ 


dr 


'o  ' dZf  'zf=i 


d~r. 


(2.23) 


For  constant  wall  temperature,  the  right-hand-side  of 


equation  (2.23)  is  one  and  R{x) -x}'^ . Hence,  any  mechanism 

\ 

to  decrease  with  time  will  yield  a 


that  causes 


growth  rate  R{z)<z 


1/2 


Substituting  the  above  dimensionless  variables  into 
the  governing  equations,  one  obtains  for  the  microlayer 


dr  Rfj  dr 


de 


n cir 


€ \ _ ^t^d L'(r)  ^ ^(^d 


Rlrdr^  dr!  R^  Z(7)  r dz^  dz}  (2.24) 


for  Osrsl,  Osz.^l.  In  the  above. 


= a/7c/  Z(7)<Z{1)  =(/>i?^(T)<i?/T). 


dr 


dt 


For  very  small  (p,  the  right-hand-side  (RHS)  of  equation 


(2.24)  reduces  to 


^d 

L\7) 


d 1-30 


+ 


dr  \ dr  I r dz^ 


^d  ^0{ 


L\?)  d-z 


and  consequently  equation  (2.24)  becomes,  to  the  leading 


21 


order  for  small  <p 


de^  _ ^ 


d^9 


dr 


R C^Fr.  'Pr'^  dzf 


for  0<r<l,  0<Zf<l, 


where  the  prime  denotes  the  time  derivative  w.r.t. 


(2.25) 
r . If 


Ci^Pr^<l,  the  conduction  term  dominates  and  a quasi-steady 
heat  conduction  problem  results.  Only  for  some  cases,  such 
as  low  Jacob  number  boiling,  it  is  found  that  Ci^Prje<l  is 
possible.  Hence  all  three  terms  on  equation  (2.25)  must  be 
retained  for  numerical  accuracy. 

It  is  also  noted  that  since  the  second  order  derivative 


wrt  r disappears  following  the  foregoing  scaling  analysis. 


only  one  boundary  condition  in  the  r'-direction  is  required. 
Because  ^^(r)>0,  the  first  order  derivative  w.r.t  r acts  as 
a convection  term  with  a negative  speed  in  the  2*-direction . 
Thus  the  condition  at  r=0  given  by  equation  (2.9)  is 
dropped.  The  advantage  of  solving  equation  (2.25)  over 
equation  (2.24)  is  that  equation  (2.25)  is  parabolic, 
instead  of  elliptic,  and  marching  from  r=l  to  r=0  for 
equation  (2.25)  is  relatively  simple. 

For  the  solid  layer,  the  energy  equation  in 
dimensionless  form  is 


dr  R^  dr  ? 8?'  dr  I fF  kp.  '^sato  (2.26) 


for  Oir<'»,  -I'S  z ^0.  For  the  first  term  on  the  RHS  of 
equation  (2.26),  the  ratio  (Rb(r)/H)^  is  small  at  the  early 
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stage  of  growth  but  can  be  large  at  a later  time  if  a thin 
heater  is  used.  Thus  unlike  the  liquid  microlayer  case,  both 
the  derivatives  wrt  r and  are  important  in  equation 

(2.26) .  Using  equations  (2.1),  (2.18),  and  (2.22),  equation 

(2.26)  can  be  expressed  as 

^ . Fo^  * 

dr  dr  a (^Ja  r dri  ^ATsato  (2.27) 


where 


Oi  t. 

s a 


(2.28) 


Fo  is  the  solid  Fourier  number.  The  last  term  in  equation 

(2.27)  is  usually  insignificant  in  most  boiling  systems.  In 
all  the  data  analyzed,  when  neglected,  the  last  term  causes 
at  most  10%  error  for  ^(r)in  a worst  case  with  large  t^. 
Although  this  term  is  kept  in  the  computation  to  obtain 
better  values  for  C^;  no  particular  attention  will  be  given 
to  its  effect  on  the  growth  rate. 

It  is  emphasized  here  that  the  introduction  of  the 
bubble  departure  time  scale,  t<j,  in  the  computation  of  the 
growth  rate  is  an  artifact.  The  growth  rate  R(t)  in 
dimensional  form  is  independent  of  t^.  However,  t^  is  a 
useful  quantity  to  characterize  the  effective  thickness  of 
the  heater  in  defining  Fo  and  to  scale  the  growth  time.  If 
td  is  quite  small  so  that  Fo<l,  the  heater  is  considered  to 
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be  very  thick  and  the  temperature  at  the  bottom  of  the 
heater  is  not  affected  by  the  time  the  bubble  departs.  If 
Fo-O(l)  or  larger  the  temperature  on  the  bottom  of  the 
heater  will  be  influenced  by  the  bubble  growth  on  a time 
scale  of  td«  The  thermal  interaction  between  the  bubble  and 
the  heater,  which  determines  the  bubble  growth  rate,  thus 
depends  strongly  on  Fo.  From  another  perspective,  the 
actual  physical  process  does  depend  on  t^  because  after 
departure  the  solid  temperature  field  will  rise  again.  An 
estimate  for  t^  can  be  obtained,  for  example,  from  the 
balance  of  the  buoyancy  force  with  the  growth  force  as 
discussed  in  Zeng  et  al.  (1993a).  In  estimating  t<j  and 
calculating  the  growth  force,  the  bubble  growth  rate,  which 
is  the  solution  to  be  sought  here,  is  needed.  This  dilemma 
can  be  easily  solved  since  we  do  not  need  the  exact  value  of 
td  in  equations  (2.25)  and  (2.27),  and  the  growth  rate  can 
be  estimated  by  using  the  constant  wall  temperature  solution 
for  R(t)  given  by  equation  (2.18).  Hence  the  estimate  used 
here  for  t<j  is 


which  gives 


Ja 


(2.29) 


^2/3  /4/3(Q  • (2.30) 


It  is  noted  that  the  above  choice  for  t^  is  adopted. 
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Non-dimensionalization  of  the  boundary  condition 
(2.14b)  on  the  liquid-solid  interface  results  in 

_ K u?) 

dz^  H dz^ 

Using  L{r)=<pr,  r=CR{t)r=R^{t^)  C equation  (2.21)  for 

^c\  ^d) 

<p  and  equation  (2.18)  for  R<,(t),  the  above  can  be  expressed 
as 


dQf  /yl/2  _ 

-4-=  — Pr\'  R(t)  at  z=0 

dz^  K dz^ 


(2.31) 


where  the  liquid-to-solid  conductivity  ratio  and  thermal 


diffusivity  ratio. 


K = kjk^  (2.32) 

a = — 

and  oLs  (2.33) 

are  two  more  independent  dimensionless  parameters  that  are 
important  to  the  bubble  growth  process  in  addition  to  Ja, 

Pr^,  and  Fo.  It  is  noted  that  k only  appears  in  the 
boundary  condition,  equation  (2.31),  while  a appears  in  both 
the  boundary  condition  and  the  energy  equation  for  the 
solid.  Hence,  the  effect  of  a on  the  growth  rate  and  the 
solid  thermal  field  is  more  complicated  than  that  of  k . 

2. 2. 2. 2 Coordinate  transformations 

To  apply  the  boundary  condition  at  large  r for  the 
solid  temperature,  the  following  simple  transformations  are 
applied  in  the  P-direction 
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r = r {l-5'tan-^[(l-C)tan(l/5')]}  forO<^<\, 

(2.34) 

where  r„  is  the  radial  position  at  which  equation  (2.12b)  is 
applied  and  the  constant  Sj.  determines  the  percentage  of  the 
grid  number  allocated  to  the  region  of  i^l . Typically, 

Sr=0.7  to  0.8  is  used  and  ranges  from  5 to  25.  To  place 
more  grids  near  z^=0  where  a large  temperature  gradient 
occurs,  the  following  transformation  is  applied 

z = -1  + 5"  tan’^[(l-e)tan(l/5'^]  forQ<s<\ 

' (2.35) 

where  is  a constant  characterizing  the  stretching.  The 

clustering  of  the  grids  occurs  near  . For  a thicker 

plate,  82=0.65  is  typically  used.  For  a thin  solid  layer,  a 

larger  value  of  82=1.0  is  preferred  because  the  thermal 

layer  will  quickly  penetrate  to  the  bottom  soon  after 

incipience . 

The  computation  must  start  from  some  very  small  but 
nonzero  initial  time,  Xq,  with  the  present  formulation. 
Because  the  error  for  R{x)  associated  with  Zq  is  of 
0(  (Zq)^^^)  , it  is  desirable  to  use  Tq  as  small  as  possible.  To 
this  end,  the  following  transformation  is  introduced 

T = 

(2.36) 

with  a uniform  "time  step"  Aa. 

The  energy  equation  for  the  solid  heater  is  transformed 
to  (^,£)  coordinates.  The  equation  is 
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1 ac  _ 1 1 a ac  \ dc 

2a  da  ^(j)  ^ dr^  C^ aJa^ 

I FO^  ( ^ \ ^ 

^de\d€  dl^  I dl^  ir>  ■ 

It  is  solved  using  the  ADI  method  (Anderson  et  al.,  1984) 
while  that  for  the  microlayer  is  solved  by  marching  from 
r=l  to  r=0  . 

2. 2. 2. 3 Solution  procedure 

Due  to  the  coupling  between  the  solid  and  liquid 
temperature  fields  through  the  boundary  conditions,  9g  and 
Gji  must  be  solved  simultaneously  at  a given  instant. 

However,  6^  needs  to  be  solved  only  for  1 while  0g  should 
be  solved  for  r'-r'„  with  . Hence,  an  iterative  solution 

procedure  is  adopted  which  is  precisely  described. 

The  solution  for  9^  and  9^  are  assumed  to  be  known  at 
a=a"  as  and  e/.  At  a=a"^^=a"+Aa,  the  first  iteration  for 
9g^*^  is  solved  using  the  ADI  scheme  by  simply  marching  in 
time  by  one  time  step  for  the  given  boundary  condition 
9s(z^  = 0)  . The  liquid  temperature  9ji  is  solved  next  with 
0"(i^=O).  After  one  sweep,  the  interface  temperature 
e"^^(i^=0)  is  updated  by  forcing  the  continuity  of  the  heat 
flux, 

dz  kji  dz 

The  second  iteration  for  begins  with  this  updated 

condition  at  z = 0 . Several  iterations  bring  the  successive 

s 
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difference  in  temperature  on  the  interface  down  to  10"^ 

(this  tolerance  is  sufficient  for  the  present  computation). 
After  the  temperature  fields  converge  at  (or  the 

growth  rate  ^(r"^^)is  updated  by  integrating  equation  (2.23) 
using  Euler's  explicit  scheme.  The  information  for 
is  a necessary  input  in  equations  (2.25)  and  (2.27). 

Although  the  solution  for  ^(r"^^)  is  only  first  order 
accurate  in  La,  the  0(Aa)  accuracy  is  not  a concern  here 
because  a very  small  La  has  to  be  used  to  retain  the 
numerical  stability  in  solving  equations  (2.23),  (2.25)  and 
(2.27).  Typically,  n=10^  is  needed  for  rsl  to  achieve  the 
convergence  of  the  solution  with  a moderate  and  large  Ja  if 
16  to  26  grids  are  allocated  in  O^rSl.  For  Ja-1,  the  total 
number  of  time  steps  is  on  the  order  of  n=8xl0^  for  Fo=10, 
100  and  1000. 


2 . 3 Results  and  Discussions 

2.3.1  Determination  of  the  Empirical  Constants  C and  C, 

After  examining  a limited  number  of  photographs  for 
growing  vapor  bubbles  in  a saturated  pool  boiling  (Van 
Stralen  et  al.,  1975,  Akiyama,  1969,  and  Hospeti  and  Mesler, 
1969),  the  following  trend  is  observed; 

(a)  For  high  Jacob  number,  the  bubble  shape  tends  to 
resemble  a hemisphere  on  top  of  a wedge  (microlayer).  Thus  C 
is  close  to  1 . 

(b)  For  low  Jacob  number,  the  bubble  shape  is  more 
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spherical  and  C is  thus  small. 

The  values  of  C are  determined  from  the  vapor  bubble 
profiles  and  the  photographs  provided  by  Van  Stralen  & Cole 
(1975),  Akiyama  et  al.  (1969)  and  Hospeti  and  Mesler  (1969). 
The  following  correlation  for  C best  fits  these  data, 

C=  [ (0.4134 /a ® + (1-0.1 

(2.38) 

Fig.  2.2  compares  the  above  fit  with  the  data.  It  is 
possible  that  C may  also  depend  on  other  parameters  such  as 
Pr^  and  Fo.  For  the  purpose  of  developing  a model  for 
predicting  the  bubble  growth  rate,  a simple  form  for  C is 
adopted  and  is  allowed  to  vary  with  Ja,  Pr^  and  a.  As 
will  be  demonstrated,  this  approach  turns  out  to  be  quite 
useful . 

The  determination  of  C^,  which  describes  the  microlayer 
wedge  angle  as  shown  by  equation  (2.5),  is  crucial  to  the 
success  of  the  present  study.  Because  the  hydrodynamics  are 
not  considered  in  this  analysis,  must  be  determined 
empirically.  For  each  experimentally  measured  R(t),  a series 
of  values  for  Cj  are  used  in  the  numerical  computation  for 
0g,  0^  and  R(t).  The  one  that  results  in  the  closest 
agreement  between  the  measured  and  computed  R(t)  over  a 
large  range  of  t prior  to  departure  is  taken  to  be  the  value 
of  Cj  required  in  the  present  analysis.  The  results  of  such 
are  given  in  Table  2.1.  The  data  are  presented  from  large 
to  small  Ja.  Based  on  these  32  data  sets,  the  following 
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Ja 


Figure  2.2  Dependence  of  C on  Ja  based  on  the  data  of  Van 
Stralen  et  al.  (1975),  Akiyama  et  al.  (1969),  and  Hospeti  and 
Mesler  (1969) . 


Table  2.1  Microlayer  angle  parameter  c.  determined  from  various  sources 
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correlation  for  is  obtained 


C -0.00643 /a 


-0.390  0.0746^-0.221 

7 lx  (Jl 


(2.39) 


by  requiring  the  least  relative  deviation 


(2.40) 


in  which  Ci*„  is  the  value  determined  from  matching  the 
experimental  growth  rate  data  with  the  present  computation. 
The  above  expression  gives  r.d.=0.232.  A very  similar 
expression  is  obtained  with  a standard  deviation  of  0.19  for 
the  exponent  of  Pr^  using  a least  squares  fitting.  Hence, 
within  statistical  error,  the  exponent  for  Pr^  may  be  taken 
as  -0.5.  Requiring  the  least  relative  deviation  again. 


with  a relative  deviation  of  0.238.  Fig.  2.3  compares  the 
above  correlation  for  with  the  data  presented  in  Table 
2.1;  reasonable  agreement  is  observed.  An  important 
implication  immediately  follows  when  one  substitutes  (2.41) 
into  (2.21),  (2.18),  (2.25),  and  (2 .27 )-(2 .31 ) , which  are 
the  only  places  where  Pr^  appears.  Since  Cj  and  Pr^  appear  in 
the  form  of  C^Pr^^^^,  it  is  obvious  that  Pr^  is  no  longer 
relevant  to  the  present  problem.  This  implies  that  the 
liquid  viscosity  is  not  important  in  determining  the  bubble 
shape  which  are  characterized  by  C and  in  the  present 
analysis . 


(2.41) 
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Figure  2 . 3 Comparison  for  based  on  the  correlation  and 
data  presented  in  Table  2.1. 
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To  understand  the  independence  of  vapor  bubble  shape  on 
the  liquid  viscosity,  it  is  instructive  to  estimate  the 
Reynolds  number  of  the  liquid  motion  associated  with  the 
bubble  growth.  Using  d=2R(t)  as  the  length  scale,  ^(t)  as 
the  velocity  scale,  assuming  R(t)-J3t^^^ , the  Reynolds  number 
of  the  liquid  flow  may  be  estimated  as  Re^~’P^/vn.  For  case 
No.  3 (Ja=766.2)  listed  in  Table  2.1,  using  the  measured 
growth  rate  shown  in  Fig.  2.4b,  it  is  easily  shown  that 
Red-12,941  for  /3^-8.8  (mm^/ms)  and  v^=6.8xl0‘^  (cm^/s). 

Hence,  the  liquid  motion  induced  by  the  expanding  bubble  in 
this  case  is  clearly  of  inviscid  nature.  Even  for  case  No. 
28  (Ja=11.24)  for  which  )8^-0.064  (mm^/ms),  one  obtains 
Red-256  by  taking  v^=5xl0’^  (cm^/s).  Typically,  Red-40  may  be 
the  lower  limit  above  which  a uniform  liquid  flow  around  a 
bubble  can  be  treated  as  inviscid  flow.  For  the  present 
case,  the  bulk  liquid  motion  is  mainly  in  the  radial 
direction.  Because  the  viscous  effect  associated  with  the 
radial  motion  is  weaker  than  that  in  a uniform  flow  over  a 
sphere,  the  liquid  flow  induced  by  the  growing  bubble  in  the 
present  study  is  mainly  of  inviscid  nature.  Therefore  the 
hydrodynamics,  which  is  important  to  the  bubble  shape 
appears  to  have  a universal  effect  on  the  bubble  growth 
rate,  i.e.  independent  of  the  liquid  viscosity,  for  most  of 
the  data  listed  in  Table  2.1.  The  foregoing  discussion 
also,  to  certain  extent,  justifies  the  choice  of  using  Ja 
alone  to  correlate  c as  given  by  equation  (2.38). 
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2.3.2  Comparison  of  the  Predicted  and  Measured  Growth  Rate 
To  validate  the  present  numerical  analysis  and  also  to 
examine  the  effect  of  the  Jacob  number,  Ja,  on  the  growth 
rate,  comparisons  are  made  between  the  present  prediction 
and  the  experimentally  measured  R(t)  over  a large  range  of 
Ja.  Fig.  2.4a  shows  such  a comparison  in  the  high  Jacob 
number  range,  Ja=201~1973 . 6 at  k=0 . 00156-0 . 0016 . The  data 
are  taken  from  Van  Stralen  et  al.  (1975)  who  used  water  as 
the  boiling  liquid.  The  Prandtl  number,  Pr^,  in  Fig.  4a 
changes  from  2.97  to  4.65;  but  it  is  not  expected  that  Pr^ 
will  affect  the  growth  rate.  Excellent  agreement  is  observed 
between  the  prediction  and  the  measurement  in  this  large  Ja 
range.  Fig.  2.4b  compares  the  numerical  results  with  the 
data  of  Cole  & Shulman  (1966)  in  which  Ja  ranges  from  90.6 
to  766.2,  k=0.027  to  0.0297  and  Pr^  =2.17-4.54.  For  these 
intermediate  to  high  values  of  Ja,  excellent  agreement  for 
R(t)  is  again  observed.  Fig.  2.4c  compares  the  prediction 
with  the  data  obtained  by  Cole  & Shulman  (1966)  using 
methanol  as  the  boiling  liquid  which  has  a higher  Prandtl 
number,  Pr^  =4.6-6.13.  The  Jacob  number  ranges  from  37.6  to 
130.0  which  is  considered  to  be  in  the  intermediate  range. 
Very  good  agreement  is  seen  between  these  data  and  the 
prediction.  While  the  measurements  shown  in  Figs.  4a-c  were 
obtained  under  subatmospheric  or  atmospheric  conditions. 

Fig.  2.4d  compares  the  present  prediction  with  the  data  of 
Staniszewski  (1959)  obtained  at  pggt=193  and  276  kpa  (or  28 
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Figure  2.4a  Comparison  of  the  bubble  growth  rate  R(t)  between 
the  present  numerical  prediction  and  the  data  of  Van  Stralen 
et  al.  (1975)  using  water  as  boiling  liquid. 
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Figure  2.4b  Comparison  of  the  bubble  growth  rate  R(t)  between 
the  present  numerical  prediction  and  the  data  of  Cole  and 
Shulman  (1966)  using  water  as  boiling  liquid. 
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Figure  2.4c  Compariison  of  the  bubble  growth  rate  R(t) 
between  the  present  numerical  prediction  and  the  data  of  Cole 
and  Shulman  (1966)  using  methanol  as  boiling  liquid. 
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Figure  2.4d  Comparison  of  the  bubble  growth  rate  R(t)  between 
the  numerical  prediction  and  the  data  of  Staniszewski  (1959) 
using  water  as  boiling  liquid  at  elevated  pressure. 
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and  40  psi),  with  Ja=11.24  and  9.16.  Even  at  this  elevated 
pressure  with  relatively  low  Jacob  number,  good  agreement 
can  still  be  observed  for  most  of  the  time.  From  Figs. 
2.4a-d,  it  is  very  obvious  that  R(t)  strongly  depends  on  Ja; 
an  increasing  Ja  results  in  an  increasing  growth  rate — a 
well  established  fact. 

It  is  also  noted  that  at  276  kpa  pressure,  the  present 
model  over-predicts  R(t)  at  the  later  stage  of  growth.  The 
data  (No.  31  and  No.  32  in  Table  2.1)  by  Akiyama  et  al. 
(1969)  were  obtained  at  even  higher  pressure.  For  those  two 
cases  the  Reynolds  number  is  of  0(1),  especially  at  the 
later  stage  of  growth.  Hence,  detailed  hydrodynamics  may 
need  to  be  properly  accounted  for  in  order  to  improve  the 
prediction.  Another  source  of  error  for  R(t)  at  low  Jacob 
number  may  be  associated  with  equation  (2.4)  in  which  the 
energy  transfer  at  the  vapor  dome  between  the  bubble  and  the 
bulk  liquid  is  neglected  based  on  the  assumption  that  bulk 
liquid  is  at  saturated  conditions.  This  would  only  cause  a 
small  error  if  R(t)  is  much  larger  than  the  thermal  layer 
thickness,  5thf  the  bulk  liquid  for  most  of  the  growth 
period.  At  low  Ja,  the  bubble  is  typically  small  and, 
hence,  the  bubble  growth  may  occur  largely  within  the 
thermal  layer.  Under  such  conditions,  the  unsteady  heat 
transfer  between  the  vapor  dome  and  the  bulk  liquid  is  no 
longer  negligible.  From  this  view  point,  the  present 
analysis  may  not  be  directly  applied  to  very  low  Jacob 
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number  boiling  systems  without  further  study. 

2 . 4 Summary 

A numerical  analysis  is  formulated  to  study  the  growth 
of  vapor  bubbles  in  saturated  heterogeneous  boiling.  The 
thermal  interaction  of  the  temperature  fields  in  the  thin 
liquid  microlayer  and  in  the  solid  heater  must  be  properly 
accounted  for  in  order  to  accurately  predict  the  growth  rate 
of  vapor  bubbles  attached  to  a heating  surface.  The  vapor 
bubble  shape  parameter,  C,  and  the  microlayer  wedge  angle 
parameter,  C^,  are  two  empirically  determined  constants 
which  are  used  in  the  numerical  computation.  In  particular, 
the  wedge  angle  parameter,  C^,  is  obtained  by  matching  the 
measured  and  predicted  bubble  growth  rate.  It  is  important 
to  note  that  once  Cj  is  correctly  specified  in  the 
computation,  the  predicted  bubble  growth  rate  compares  very 
well  with  the  reported  experimental  data  over  a wide  range 
of  conditions  during  the  entire  growth  period.  The 
empirically  determined  correlation  for  Cj  suggests  that  the 
liquid  Prandtl  number,  hence  the  viscosity  of  the  liquid, 
plays  little  role  in  the  bubble  growth  process.  A 
subsequent  estimate  indicates  that  the  Reynolds  number  of 
the  bulk  liquid  motion  induced  by  the  bubble  is  of  inviscid 
nature;  thus  the  viscosity  is  not  important  in  determining 
the  parameters  C and  Cj. 

The  dimensional  analysis  for  the  governing  equations 
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and  boundary  conditions  shows  that  there  exists  four 
relevant  dimensionless  parameters  for  bubble  growth:  Jacob 
number  Ja,  Fourier  number  Fo  of  the  solid  based  on  the 
thickness  of  the  heater  and  the  bubble  departure  time  scale, 
the  liquid-to— solid  conductivity  ratio  K , and  the 
liquid-to-solid  thermal  diffusivity  ratio  a.  The  Jacob 
number  is  the  most  important  one  affecting  C,  C^,  and 
dimensionless  growth  rate  J?(t)  . The  effects  of  the  four 
parameters  on  J?(z)  and  the  thermal  fields  of  the  liquid  and 
the  slid  heater  will  be  discussed  in  chapter  3. 


CHAPTER  3 

VAPOR  BUBBLE  GROWTH  IN  SATURATED  POOL  BOILING 
— PART  II.  GROWTH  RATE  AND  THERMAL  FIELDS 

3 . 1 Background 

In  Chapter  2 of  the  present  study,  the  governing 
equations,  the  boundary  conditions,  and  the  initial 
conditions  for  the  unsteady  energy  transfer  from  the  solid 
heater  through  the  liquid  microlayer  to  the  vapor  bubble 
have  been  systematically  formulated.  A time  varying 
coordinate  system  is  used  to  resolve  the  large  temperature 
variations  over  a small  area  beneath  the  vapor  bubble.  The 
bubble  shape  parameter,  C,  was  correlated  with  Jacob  number 
based  on  experimental  bubble  shape  profiles,  and  the 
microlayer  wedge  angle  parameter,  C^,  was  determined  by 
matching  the  finite  difference  solution  for  the  growth  rate, 
R(t) , with  the  existing  experimental  data  over  a large  range 
of  conditions.  The  numerical  results  for  the  bubble  growth 
rate  have  compared  very  well  with  the  existing  data  for  the 
growth  rate  over  a wide  range  of  conditions.  From  the  basic 
formulation,  qualitative  analysis,  and  the  result  for  C^ , it 
is  observed  that  the  Prandtl  number  has  little  effect  on  the 
growth  process. 

The  problem  of  vapor  bubble  growth  in  saturated 
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heterogeneous  boiling  involves  many  physical  properties  and 
boiling  conditions:  h^g,  p^,  p^,  k^,  k^,  a^,  a^,  p^, 

(or  , and  H (see  nomenclature)  which  appear  in  the 

governing  equations  for  the  three  phases,  boundary 
conditions,  and  initial  conditions.  The  dependent  variables 
are  R(t)  , T^,  and  T^.  The  liquid  kinematic  viscosity  is 
listed  above  formally  because  in  general  the  bubble  shape 
depends  on  the  hydrodynamics;  it  appeared  in  a previous 
investigation  (Cooper  and  Lloyd,  1969) . For  intermediate 
and  high  Jacob  numbers,  it  has  been  shown  in  Chapter  2, 
based  on  a qualitative  estimate  for  the  Reynolds  number  and 
an  examination  of  existing  bubble  growth  rate  data,  that  p^^ 
is  not  important.  No  particular  attention  was  paid  to  the 
explicit  effects  of  q^^  (or  q^^^)  on  R(t)  because  the 
contributions  are  usually  small  during  the  majority  of  the 
growth  period.  However,  the  implicit  effects  of  q^^  (or  q^^^) 
through  AT^^^q  R(t)  are  important  and  were  incorporated. 

Further  examination  of  the  dimensionless  forms  of  the 
governing  equations,  boundary  and  initial  conditions 
indicated  that  the  dimensionless  growth  rate 

mr)  = ^(0  . ° 

■^.(’■=1)  (3.1) 

only  depends  on  four  dimensionless  parameters:  Jacob  number 
Ja,  Fourier  number  Fo,  liquid-to-solid  conductivity  ratio  k, 
and  liquid-to-solid  diffusivity  ratio  a.  In  the  above,  R^ 
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is  the  solution  for  the  growth  rate  based  on  a constant  wall 
temperature  and  t^j  is  the  time  scale  for  bubble  departure 
estimated  from  the  analysis  of  Zeng  et  al.  (1993a).  Hence 
the  complexity  in  understanding  and  describing  the  bubble 
growth  rate  is  significantly  reduced. 

In  order  to  further  understand  the  dynamics  governing 
vapor  bubble  growth  in  heterogeneous  boiling,  a systematic 
investigation  is  presented  which  considers  the  effects  of 
(Ja,  Fo,  K,  a)  on  the  growth  rate  R{z)  and  the  thermal 
fields  of  the  liquid  and  solid.  The  numerical  results  show 
that  the  liquid  temperature  in  the  microlayer  is  practically 
linear  for  a wide  rage  of  conditions.  The  analysis  and 
numerical  solution  have  elucidated,  for  the  first  time,  the 
development  of  the  unsteady  thermal  field  in  the  heater 
during  a very  short  period  of  time  through  a very  small 
volume  which  resides  beneath  a growing  vapor  bubble. 
Experimental  techniques  described  in  the  existing  literature 
are  not  capable  of  resolving  the  fine  temporal  and  spatial 
scales  considered  investigation. 

3 . 2 Effect  of  Ja,  Fo,  k,  and  a on  R(z) 
since  experimental  measurements  which  reveal  the  effect 
of  K (or  a,  Fo)  while  holding  Ja  and  other  parameters 
constant  are  not  available,  the  parametric  study  is 
especially  useful  for  further  understanding  of  the  bubble 
growth  process.  In  what  follows,  the  effects  of  varying  Ja, 
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Fo,  K,  and  a will  be  examined  separately. 

Fig.  3.1a  shows  the  effect  of  varying  Ja  (from  1 to 
1000)  on  i?(r)  at  Fo=l,  k=0.005  and  a=0.005.  Because  RcCt^) 
strongly  depends  on  Ja  as  indicated  by  equations  (2.18)  and 
(2.38),  the  most  significant  part  of  the  Jacob-number-ef feet 
on  R(t)  has  been  accounted  for  by  R^,(tjj)  . Therefore  the 
difference  in  R (v)  among  different  values  of  Ja  in  Fig.  3.1a 
mainly  reflects  the  effect  of  the  coupling  of  the  thermal 
fields  between  the  solid  heater  and  liquid  microlayer.  It 
will  be  shown  in  section  3.3  that  the  liquid  temperature 
profile  in  the  microlayer  is  practically  linearly.  Hence 
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dz 


practically  does  not  vary  with  during  bubble  growth. 


^ 1 yO  n *7  R 0 

On  the  liquid-solid  interface,  since  C^Pi^  <x.  Ja'  , a small 
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Ja  implies  a small  — for  a given  , which  results  in  a 


3z. 


3z 


smaller  R(r),  as  seen  from  equations  (2.23)  and  (2.31). 

Because  R (r)  also  appears  in  equation  (2.31),  a smaller  R (r) 

39^ 

is  consistent  with  a smaller  • 

Fig.  3.1b  shows  the  effect  of  changing  Fo  (from  0.01  to 
10000)  for  Ja=10,  k=0.005,  and  a=0.005  for  T=t/t^<l.  Little 
difference  in  R (r)  between  Fo=0.01  (a  thicker  heater)  and 
Fo=l  is  observed  for  this  set  of  parameters.  As  Fo 
increases,  say  due  to  decreasing  thickness  of  the  heater,  H, 
or  longer  departure  time,  t^,  R (r)  decreases.  The  decrease 
in  H implies  a reduction  of  the  heater  thermal  capacity 
which  will  reduce  the  growth  rate.  Under  otherwise 
identical  conditions,  a larger  value  of  t^  implies  more 
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Figure  3 . la  Parametric  dependence  of  the  normalized  growth 
rate,  Ja=l,  10,  100  and  1000  at  (Fo, /c, a)  = ( 1 , 0 . 005 , 0 . 005)  . 
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Figure  3 . lb  Parametric  dependence  of  the  normalized  growth 
rate,  Fo=0.01,  1,  10,  100,  1000  and  10000  at  ( Ja , k, a)  = ( 10 , 
0.005,  0.005) . 
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energy  is  taken  from  the  heater  between  t=0  (t=0)  and  t=t^j 

(t=1) . Hence  the  growth  rate  at  later  stages,  when  t is 

close  to  1,  is  slower  for  larger  t^^  in  comparison  with  that 

for  smaller  t..  To  better  understand  the  effect  of  Fo  on 

the  growth  rate,  consider  two  cases:  1)  Fo^=l  and  2)  Fo2=10 

while  all  other  parameters  are  held  the  same.  If  the 

computation  for  stops  at  t^=10  while  that  for  R2  stops  at 

Ty  =1,  the  results  for  the  two  dimensional  growth  rates  are 

identical  if  the  same  dimensional  time  is  used,  i.e., 

Fo, 

R^  (Fo^t,)=R2(Fo2T2)  if  This  simply  means  that  if 

FOi 

one  is  interested  in  R(t)  for  0<t<<»  the  use  of  Fo  is 

unnecessary  and  the  relevant  time  scale  should  be 

a^t 

Instead,  an  instantaneous  Fourier  number  Fo  = —^  should  be 

H ^ 

used  to  characterize  the  thermal  interaction.  In  one-g 
boiling,  the  bubble  departs  at  some  finite  time  t^;  thus  no 
experimental  information  is  available  for  R(t)  beyond  t^j. 
Since  it  is  always  desirable  to  sufficiently  resolve  the 
growth  rate  between  t=0  and  t=t^,  the  use  of  T=t/t^  and  Fo 
makes  the  computation  and  the  interpretation  of  R(t) 

dR 

considerably  easier.  It  is  also  noted  that  as  t-«»,  — -*  0 

dt 

because  there  is  insufficient  energy  in  the  heater 
underneath  the  microlayer  to  sustain  the  growth.  This 
scenario  is  equivalent  to  the  case  when  Fo-«o,  as  was 
observed  in  the  computations. 

Fig.  3.1c  shows  the  effect  on  R (r)  of  varying 
conductivity  ratio  k (from  0.0005  to  0.05)  at  Ja=10,  Fo=l, 
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Figure  3.1c  Parametric  dependence  of  the  normalized  growth 
rate,  /c=0.0005,  0.001,  0.01  and  0.05  at  (Ja,  Fo,  a)  = (10,  1, 
0.005) . 
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and  a=0.005.  Since  k appears  only  in  the  boundary  condition 
(2.31),  it  is  easy  to  see  that  a larger  k (which  may  result 

from  using  a low-conductivity  solid  heater)  implies  a 

60.  60,  60, 

smaller  — for  a given  — . Of  course,  is  not  fixed; 

dz^  6z,  ^ 6z,  ^ ^ 

it  increases  with  increasing  k as  will  be  shown  in  Fig.  3.7. 

1 30, 

Nevertheless  the  net  effect  of  increasing  k such  that -::r 

de,  . . 

decreases.  A smaller  naturally  implies  a smaller  amount 

of  energy  is  conducted  to  the  bubble  which  results  in  a 

smaller  growth  rate.  Hence  -R  (r)  decreases  as  k increases 

and  vice  versa,  as  shown  in  Fig.  3.1c. 

Fig.  3.  Id  shows  the  effect  on  R(z)  of  changing  a 

(0.0005  to  0.05)  at  Ja=10,  Fo=l,  and  k=0.005.  Contrary  to 

the  effect  of  k,  an  increasing  a results  in  an  increasing 

30, 

R (t)  . However,  equation  (2.31)  indicates  that  — , hence 
- 

R(r),  increases  with  increasing  From  Figs.  3.1c-d,  it 

is  seen  that  when  k or  a changes,  say,  by  a factor  of  10, 
the  resulting  change  in  R (r)  is  larger  in  Fig.  3.1c  than  in 
Fig.  3.  Id.  This  simply  results  from  the  factor  in 

equation  (2.30). 

To  recapitulate,  the  effects  of  Ja,  Fo,  k and  a on  the 
normalized  growth  rate  R (r)  are:  1)  increasing  Ja  and  a will 
result  in  an  increasing  R(t);  2)  increasing  Fo  and  k will 
result  in  a decreasing  R(r). 


3 . 3 Thermal  Field  of  the  Liquid  Microlaver 
It  is  not  clear  from  equation  (2.25)  what  the  shape  of 
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Figure  3 . Id  Parametric  dependence  of  the  normalized  growth 
rate,  a=0.0005,  0.001,  0.001  and  0.05  at  (Ja,  Fo,  k)=(10,  1, 
0.005) . 
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the  liquid  microlayer  temperature  6.  will  be  in  the 


z-direction  under  general  conditions.  Only  when  C^Pr^  is 

small,  it  is  clear  that  the  temperature  will  be  almost 

linear  in  ^ because  diffusion  the  term  in  (2.25)  is  clearly 

dominant.  In  order  to  qualitatively  assess  the  linearity  of 

6f  from  a linear  profile,  we  calculate  the  change  in  the 

d9, 

slope,  relative  to  the  slop.  Since  the  change  in  the 

slope  can  be  conveniently  characterized  by  the  second  order 

d^e. 


derivative. 


dzl 


, we  use  the  following  parameter 


rw  = i E 

^ ’ n=l 


dZf  ' 


(3.2) 


to  represent  the  deviation  from  the  linear  profile.  For 
r«l  the  deviation  of  0^  from  a linear  profile  is  quite 
small  and  as  F-^0,  the  profile  is  exactly  linear.  In  the 
above,  the  summation  E is  carried  over  N interior  grid 
points  of  the  computational  domain  within  the  microlayer. 

Fig.  3.2a  shows  the  variation  of  F during  the  bubble 
growth  period,  0<t<1,  for:  i)  (Fo,/c,  a)  = (l,  0.005,  0.005) 
with  Ja=10  and  1000;  and  ii)  (Ja,  k,  a)=(10,  0.005,  0.005) 
with  Fo=10  and  1000.  Fig.  3.2b  shows  the  dependence  of  F on 
K and  a:  i)  (Ja,  Fo,  a)  = (10,  1,0.005)  with  /c=0.05  and 
0.0005;  ii)  (Ja,  Fo,  k =(10,  1,  0.005)  with  a=0.05  and 
0.0005.  The  largest  value  of  F always  occurs  at  the  initial 
time  step.  Under  otherwise  identical  conditions,  a larger 
Ja  results  in  a larger  F.  This  is  due  to  the  fact  that 
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Figure  3.2a  Deviation  from  the  linearity  of  the  liquid 
temperature  profile  in  the  microlayer,  (Fo,/c,a)  = (l,  0.005, 

0.005)  with  Ja=10  and  1000  and  ( Ja, k, a) = ( 10 , 0.005,  0.005) 

with  Fo=10  and  1000. 
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Figure  3.2b  Deviation  from  the  linearity  of  the  liquid 
temperature  profile  in  the  microlayer,  ( Ja, Fo, a) = (10 , 1, 

0.005)  with  k=0.0005  and  0.05  and  (Ja, Fo,k) =(10,  1,  0.005) 

with  a=0.0005  and  0.05. 
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larger  Ja  results  in  a larger  wedge  angle,  or  thermal 
capacity,  of  the  microlayer,  as  shown  by  equations  (2.5)  and 
(2.38).  For  Ja=10,  F<0.01  for  t>0.001.  Even  for  Ja=1000,  F 
is  less  than  2%  for  t>0.01  (or  for  99%  of  the  growth 
period) . No  substantial  difference  in  F is  observed  as  Fo 
is  changed  from  10  to  1000;  the  values  of  F are  all  less 
than  1%.  A larger  F occurs  at  /c=0.05  (a  large  value)  because 
a large  k is  associated  with  a rapid  drop  of  the  solid 
temperature  on  the  liquid-solid  interface  underneath  the 
microlayer,  as  will  be  demonstrated  in  Fig.  3.7.  This 
larger  transient  variation  of  the  temperature  on  the 
boundary  causes  the  microlayer  temperature  to  deviate  more 
from  a linear  profile  because  the  thermal  diffusion  term  in 
equation  (2.25)  is  not  as  dominant  as  it  is  with  a smaller 
transient  variation.  Little  variation  of  F is  observed  when 
a is  changed  from  0.0005  to  0.05.  Even  for  k=0.05,  the 
largest  F is  only  about  3-4%  at  t=0.001  as  shown  in  Fig. 
3.2b.  It  also  needs  to  be  commented  that  F is  sensitive  to 
the  initial  time,  Tq,  when  the  computation  is  initiated.  In 
Fig.  3.2,  Tq=0.0005  was  used  and  the  results  are  presented 
for  r >0.001.  When  the  initial  time  Tq  is  reduced  to 
0.00001,  F at  T=0.001  is  further  reduced  by  a factor  of  2 or 
more.  The  reason  is  that  the  initial  conditions  for  the 
liquid  and  solid  temperatures  given  by  (2.16-17)  are 
inconsistent  with  the  actual  thermal  fields,  although  the 
inconsistency  does  not  affect  the  accuracy  of  the  bubble 
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growth  rate.  Hence,  the  actual  values  of  F may  be  less  than 
what  were  shown.  The  result  of  this  parametric 
investigation  is  that  the  liquid  temperature  profile  within 
the  microlayer  can  be  taken  as  linear  for  practical 
purposes.  This  result  will  be  used  as  a starting  point  in 
chapter  4 to  develop  a simplified  bubble  growth  rate  model. 

3 . 4 Thermal  Field  of  the  Solid  Heater 
To  understand  the  energy  transfer  process  in  the  solid 
heater,  a qualitative  description  of  the  solid  temperature 
field  for  a thick  heater,  Fo«l,  is  first  given  based  on  the 
governing  equations  without  resorting  to  extensive 
computations.  For  a thick  heater,  the  temperature  on  the 
bottom  of  the  solid  layer  is  unaffected  by  the  heat  transfer 
from  the  solid  to  liquid  for  t<t^.  When  energy  is  taken  away 
from  the  solid,  an  unsteady  thermal  boundary  layer  develops 
in  the  z-direction,  there  is  no  convection  term  in  (2.27). 
Hence,  the  relevant  length  scale  for  Fo«l  in  the 
z-direction  is  6 = Ja  t , which  is  a measure  of  the 
temperature  penetration  depth  in  the  z-direction.  In  the 
r-direction,  6 is  affected  strongly  by  the  movement  of  the 
bubble  (through  the  microlayer)  and  the  relevant 
"convective”  length  scale  is  thus  Over  a distance 

of  r~5.. , 6.  changes  appreciably.  It  is  noted  that  r=R.  (t) 
is  also  the  location  of  the  discontinuity  in  the  boundary 
condition  at  z =0  given  by  equations  (2.14b)  and  (2.15b.  As 
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a rough  estimate,  the  constant  wall  temperature  solution  for 
the  bubble  radius  given  by  equation  (2.18)  can  be  used  to 
characterize  S . If  dimensionless  isothermal  lines,  0 , are 

rc  » 

sketched,  the  contours  would  resemble  ellipses  with  and 
S being  the  major  and  minor  axes,  as  shown  in  Fig.  2.1. 
Since 


is  roughly  independent  of  time,  the  ^^-contours  are  expected 

to  expand  in  a form  of  similar  ellipses.  For  the  case  with 

large  Ja,  penetration  of  the  temperature  in 

z-direction  can  be  quite  shallow. 

When  H is  not  large  or  Fo~0(l) , the  foregoing 

description  on  the  qualitative  features  of  the  temperature 

contours  are  still  valid  near  the  nucleation  site  for  small 

t as  long  as  the  instantaneous  solid  Fourier  number  is 

a t 

small,  i.e.,  Fo-—^<l.  At  later  times,  however,  the 


temperature  contour  will  be  drastically  different  from  those 
sketched  in  Fig.  2.1.  For  Fo~0(l)  or  larger,  we  resort  to 
the  computational  results  in  order  to  elucidate  the  thermal 
fields  of  the  heater. 

Fig.  3.3a  shows  the  four  computed  temperature  contours 
corresponding  to  case  No. 22  in  Table  2.1  with  Ja=23.4, 
Fo=0.057,  /c=0.0073,  and  a=0. 00577,  at  four  different  times: 
t=0.05,  0.1,  0.5,  and  1.  The  parameter  C^=0.12  is  taken 


(3.3) 
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'C=0.05 


x=0.1 


Figure  3.3a  Contours  of  temperature  6 in  a thick  heater  at 
(Ja,Fo,/c,a)  = (23 .4,  0.057,0.0073,0.00577)  for  T=0.05,  0.1, 

0.5  and  1.  This  condition  corresponds  to  the  experiment  of 
Hospeti  and  Mesler  (1969).  All  dimensions  are  in  mm. 


from  Table  2.1.  These  four  times  roughly  cover  the  early 
stages,  later  stage,  and  "departure”  stage  of  the  growth 
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process.  Since  Fo<<l,  the  heater  is  quite  thick  and 
foregoing  asymptotic  description  is  valid  for  all  time.  The 
circle  marker  at  z=0  is  the  location  of  the  microlayer  edge, 
r=Ru(t) . It  is  clearly  seen  that  the  temperature  contours 
in  the  solid  penetrate  down  and  move  to  the  right  at  a 
similar  rate  as  the  bubble  grows.  Fig.  3.3b  shows  the  0^— 
contours  corresponding  to  case  No. 19  in  Table  2.1,  with 
Ja=25.3,  Fo=63.3,  /c=0. 00173,  and  a=0. 00144,  at  four 
different  times:  t=0.05,  0.1,  0.5,  and  1.  Since  Fo»l,  the 
heater  is  very  thin  and  the  asymptotic  description  is  not 
valid  for  all  times.  Nevertheless,  for  t=0.05,  the  6^- 
contours  in  a small  region  near  the  bubble  center,  r«R(j(t)  , 
still  agree  with  the  asymptotic  description.  On  a larger 
scale,  r~R[j(t)  or  r- 1 , the  0g-contours  between  the  case  of 
Fo>>l  and  Fo«l  display  a drastic  difference.  For  a very 
thick  heater,  the  heat  flux  in  the  z-direction  from  below  is 
always  maintained.  For  a thin  heater  the  thermal  capacity 
is  small.  As  t increases,  ^.-contours  penetrate  down  and 
reach  the  bottom  of  the  heater;  hence  the  thermal  field  near 
the  liquid-solid  interface,  which  determines  the  energy 
transfer  to  the  bubble  via  the  microlayer,  is  also 
influenced  by  the  bottom  of  the  heater.  The  growth  of  the 
bubble  quickly  depletes  the  energy  stored  in  the  region 
0<r<Ru(t)  and  -H<Zg<0.  To  maintain  bubble  growth,  albeit 
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x=0.05 


x=0.1 


x=0.5 


x=1.0 


Figure  3.3b  Contours  of  temperature  6^  in  a thin  heater  at 
(Ja,Fo,K,a)=(25. 3, 63. 3, 0.00173, 0.00144)  for  t=0.05,  0.1,  0.5 
and  1.  This  condition  corresponds  to  the  experiment  of 
Staniszewski  (1959) . All  dimensions  are  in  mm. 
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slow,  energy  must  be  drawn  from  the  region  outside  r=R|^(t)  . 
This  creates  an  energy  flow  in  the  negative  r-direction 
outside  r=Rb(t)  within  the  thin  heater.  The  relatively  long 
distance  the  energy  must  flow,  due  to  conduction,  certainly 
results  in  a slower  rate  of  energy  transfer  to  the  bubble 
and  hence  a slower  vapor  bubble  growth  rate.  This 
qualitatively  explains  the  discerned  decreasing  growth  rate 
with  increasing  Fo  shown  in  Fig.  3.1b.  Another  interesting 
and  useful  revelation  from  Fig.  3.3b  is  that  even  for  a very 
thin  heater,  the  temperatures  at  the  bottom  and  at  the 
interface  in  r<Rj^(t)  are  not  the  same.  Hence,  the 
temperature  measurement  on  the  bottom  of  heater  fails  to 
present  the  actual  thermal  field  near  the  liquid-solid 
interface  underneath  the  rapidly  growing  vapor  bubble. 

To  parametrically  examine  the  effects  of  Ja,  Fo,  k and 
a on  the  solid  temperature  field  during  the  growth  period, 
two  particular  locations  are  selected:  the  centerline  r^O 
and  the  liquid-solid  interface  2^=0.  Temperatures  at  four 
times  are  examined:  t=0.05,  0.1,  0.5  and  1 which  represent 
the  early  stages,  later  stage  and  "departure”  stage  of 
bubble  growth. 

3.4.1  Effect  of  Ja 

Figs.  3.4a-d  show  at  r=0 , Fo=l,  k=0.005,  and 

a=0.005  with  Ja=l,10,  100  and  1000.  The  growth  rate  for 
Ja=l  is  smaller  compared  with  that  for  Ja=1000;  hence,  the 
energy  removed  from  the  solid  is  less  for  Ja=l  than  for 
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Z/H 

(a) 


Z/H 

(b) 


Z/H 

(c) 


Figure  3.4  Effect  of  Jacob  number  Ja  on  the  solid  heater 
temperature  0^  along  the  centerline  at  various  times  (t=0.05, 
0.1,  0.5  and  1)  for  (Fo, /c,  a)  = ( 1 , 0.005,  0.005).  (a)  Ja=l;  (b) 
Ja=10;  (c)  Ja=100  and  (d)  Ja=1000. 
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Ja=1000.  While  the  temperature  near  the  bottom  is  nearly 
equal  to  one  for  Ja=l,  6^  through  the  whole  heater  (at  f^O) 
is  reduced  for  Ja=100  and  1000  due  to  a large  amount  of 
energy  removed  from  the  heater.  It  also  appears  that  the 
thermal  diffusion  in  the  z-direction  is  enhanced  with 
increasing  Ja;  although  it  is  clear  from  (2.27)  that  the 
diffusion  in  the  z-direction  is  independent  of  Ja.  To 
understand  this  phenomena  the  radial  variations  of 
6 (^=r/R  (t .) ) are  examined  next.  Figs.  3.5a-d  show  the 

S 0 o 

development  of  the  thermal  layer  in  the  radial  direction  at 
(Fo,  K,  a)=(l,  0.005,  0.005)  for  Ja=l,  10,  100  and  1000. 

Also  marked  on  each  curve  are  the  locations  of  the 
microlayer  edge,  r=R|3(t)  or  r=  1 . For  Ja=1000,  there  is 
hardly  any  thermal  diffusion  in  the  r-direction  beyond 
r=R|j(t)  during  the  entire  growth  process.  For  Ja=l,  very 
significant  thermal  diffusion  in  the  r-direction  can  be 
identified  for  r>R|^(t)  , in  contrast  to  the  very  small 
diffusion  of  6 in  the  z-direction  at  the  same  Ja  shown  in 
Fig.  3.4a. 

To  elucidate  the  effect  of  Ja  on  the  development  of 
the  thermal  field  within  the  solid,  especially  in  the 
z-direction,  equation  (2.27)  is  revisited.  The  second  term 
on  the  left-hand  side  is  the  convective  energy  transport  due 
to  the  bubble  growth.  The  first  term  on  the  RHS  is  the 
energy  transport  due  to  thermal  diffusion  in  the 
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Figure  3.5  Effect  of  Jacob  number  Ja  on  9^  on  the  liquid-solid 
interface  at  various  times  for  (Fo,ic,a)  = (l,  0.005,  0.005).  (a) 
Ja=l;  (b)  Ja=10;  (c)  Ja=100  and  (d)  Ja=1000.  Symbols  indicate 
the  edge  of  the  microlayer,  r=R|j(t)  . 
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i^(C)CiP^ 

r-direction.  For  a fixed  (Fo,  k,  a),  the  factor  , 

aC^Jc? 

hence  the  extent  of  r-dif fusion  is  large  at  small  Ja  and  it 
decreases  with  increasing  Ja.  As  Ja  increases,  the  growth 
rate  increases  so  that  more  energy  from  the  heater  is 
transferred  to  the  bubble  via  the  microlayer,  as  discussed 
previously.  Since  the  r-diffusion  is  decreasing  with 
increasing  Ja,  the  z-dif fusion  must  increase  in  order  to 
provide  the  thermal  energy  to  the  growing  bubble.  This 
explains  the  trend  in  9^(z^)  shown  in  Figs.  3.4a-c  for  Ja=l, 
10,  100.  However,  as  Ja  becomes  quite  large,  the  extent  of 
r-diffusion  is  limited.  Since  the  bubble  grows  so  fast,  the 
convective  heat  transfer  due  to  the  bubble  growth  becomes 
important  in  that  the  edge  of  the  expanding  microlayer 
always  experiences  a solid  surface  temperature  that  is  close 
to  1.  Hence,  at  high  Ja,  the  energy  transfer  to  the 
microlayer  is  accomplished  within  a relatively  thin  thermal 
boundary  layer  near  the  liquid-solid  interface;  within  this 
boundary  layer,  the  r-convection  and  the  z-dif fusion  are  the 
dominant  terms  in  equation  (2.27).  Thus,  the  diffusion  in 
the  z-direction  at  Ja=1000  does  not  penetrate  as  far  toward 
the  bottom  of  the  heater  as  in  the  case  when  Ja=100. 

3.4.2  Effect  of  Fo 

Figs.  3 . 6a-d  show  0(i)  at  r=  0 for  (Ja,  k,  a)  = (10, 

^ s s 

0.005,  0.005)  with  Fo  = 0.01,  1,  10,  and  1000.  For  Fo=0.01, 
the  bottom  portion  of  the  heater  is  unaffected  by  the  growth 
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Z/H 

(a) 


Z/H 

(c) 


Figure  3.6  Effect  of  Fourier 
centerline  at  various  times  for 
(a)  Fo=0.01;  (b)  Fo=l;  (c)  Fo=10 


Z/H 

(b) 


Z/H 

(d) 


number  Fo  on  6^  along  the 
(Ja,/c,a)  = (10,  0.005,  0.005). 
and  (d)  Fo=100. 
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of  the  bubble  during  the  entire  growth  process  while  for 
higher  values  of  Fo  the  bottom  temperature  at  r=  0 reaches 
0.77,  0.24,  and  0.033  near  departure.  The  temperature 
gradient  on  the  liquid-solid  interface  z^= 0 is  largest  for 
Fo=0.01;  the  gradient  decreases  with  increasing  Fo,  and  thus 
the  growth  rate  does  as  well.  This  is  entirely  consistent 
with  the  results  shown  in  Fig.  3.1b  for  the  effect  of  Fo  on 
the  growth  rate.  For  Fo  changing  from  0.01  to  1,  10,  100, 
1000  and  10000,  the  radial  distribution  of  (^=r/R^ (t^) ) is 
very  similar  to  that  in  Fig.  3.5b;  hence  they  are  not 
presented.  The  insensitivity  of  6^(^)  is  primarily  due  to 
the  explicit  independence  of  the  radial  diffusion  on  Fo  and 
the  proper  scaling  of  r using  . 

3.4.3  Effect  of  k 

Figs.  3.7a-d  illustrate  the  effect  of  the  conductivity 
ratio  K on  0(5)  at  r=0.  At  (Ja,  Fo,  a)  =(10,  10,  0.005),  k 

s s 

ranges  from  0.0005,  0.001,  0.01  to  0.05.  For  a very  small  k 
(say  k=0.0005  resulting  from  using  a high-conductivity 
heater) , only  a small  gradient  in  0 is  needed  to  supply 
energy  to  the  bubble  through  the  microlayer.  For  larger  k 
(say  K=0.01)  the  required  temperature  change  is  much  larger 
in  order  to  maintain  the  energy  flow  to  the  microlayer. 

Such  a description  is  consistent  with  the  results  presented 
in  Fig.  3.1c.  Figs.  3.8a-b  show  the  radial  variation  of  6^ 
at  (Ja,  Fo,  a)  = (10,  1,  0.005)  for  /c=0.001  and  0.01. 

Together  with  Fig.  3.5b  in  which  k=0.005,  it  is  seen  that  an 
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(a)  (b) 


Figure  3.7  Effect  of  the  conductivity  ratio  k on  6^  along  the 
centerline  at  various  times  for  ( Ja, Fo, a) = (10 , 1,  0.005).  (a) 
k=0.05;  (b)  K=0.01;  (c)  /c=0.001  and  (d)  /c=0.0005. 
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(a) 


(b) 


Figure  3.8  Effect  of  the  conductivity  ratio  k on  0^  on  the 
liquid-solid  interface  at  various  times  for  (Ja,Fo,a)=(10,  1, 
0.005).  (a)  K=0.001  and  (b)  k=0.01.  Symbols  indicate  the 
edge  of  the  microlayer,  r=R|j(t)  . 
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increase  in  the  temperature  underneath  the  bubble  is 
associated  with  decreasing  k,  which  is  similar  to  the  trend 
in  the  z-direction.  The  result  at  /c=0.0001  (not  presented 
here)  simply  reveals  that  0„(^=O)  is  very  close  to  0.95  for 

w 

0<T<1. 

3.4.4  Effect  of  a 

Figs.  3.9a-b  show  the  effect  of  a on  the  OJzJ.  It  is 
seen,  by  comparing  these  two  cases,  that  a small  a results 
in  a steeper  gradient  of  0^  at  z^=0  but  a lesser  penetration 
of  the  thermal  layer  through  the  thickness.  The  large 
gradient  of  0^  on  the  interface  at  small  a is  needed  to 
maintain  the  energy  flow  to  the  microlayer  during  the  growth 
of  the  bubble,  as  seen  from  equation  (2.31).  However,  a 
larger  temperature  variation  for  large  a (which  may  result 
from  a small  solid  thermal  diffusivity)  is  observed  in  the 
z-direction  along  the  centerline.  This  is  quite  similar  to 
the  effect  of  Ja  on  Figs.  3.10a-d  show  (^=r/R^ (t^j) ) 

at  (Ja,  Fo,  /c)  = (10,  1,  0.005)  for  a=0.0005,  0.001,  0.01,  and 
0.05.  An  increasing  (decreasing)  a results  in  reduced 
(enhanced)  thermal  diffusion  in  the  r-direction.  Since  the 
z-direction  conduction  is  controlled  by  Fo,  not  by  a,  the 
energy  flow  is  more  likely  to  be  in  the  radial  direction 
with  smaller  a (see  Fig.  3.10a  for  a=0.0005)  than  with 
larger  a (see  Fig.  3.10d  for  a=0.05).  Correspondingly,  the 
z-direction  diffusion  is  larger  for  a=0.05  than  for  a=0.0005 
as  shown  in  Fig.  3.9. 
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1 + z/H 

(a) 


1 + z/H 
(b) 


Figure  3.9  Effect  of  the  diffusivity  ratio  a on  0^  along  the 
centerline  at  various  times  for  (Ja,Fo,K)=(10,  1,  0.005).  (a) 
a=0.0005  and  (b)  a=0.05. 
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Figure  3.10  Effect  of  the  diffusivity  ratio  a on  6^  on  the 
liquid-solid  interface  at  various  times  for  (Ja,Fo,/c)  = (10,  1, 
0.005).  (a)  a=0.0005;  (b)  a=0.001;  (c)  a=0.01  and  (d)  a=0.05. 
Symbols  indicate  the  edge  of  the  microlayer,  r=R|j(t)  . 
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From  the  variation  of  6AzJ , it  is  also  seen  that:  a) 

s s 

there  is  a large  difference  between  0^  on  the  liquid-solid 
interface  and  6^  at  the  bottom  of  the  heater;  b)  the 
temperature  variation  on  is  typically  large;  c)  the 
variation  of  0,  depends  on  the  four  parameters  (Ja,  Fo,  k, 
a)  in  a complex  manner.  Hence,  measurements  of  the  solid 
temperature  at  the  bottom  of  the  heater  does  not  clearly 
reveal  the  thermal  field  at  the  solid-liquid  interface  just 
beneath  the  bubble. 


3 . 5 Summary 

The  numerical  results  show  that  the  liquid  temperature 
within  the  microlayer  is  practically  linear  during  the 
entire  growth  period  over  a large  range  of  conditions.  The 
normalized  growth  rate  R{z)=R(v)/R^{z-l)  clearly  shows  the 
effects  of  (Ja,  Fo,  k,  a) : increasing  Ja  or  a results  in  an 
increased  R (r ) while  increasing  Fo  or  k results  in  a 
decreasing  R{z).  The  asymptotic  description  and  numerical 
solution  have  elucidated  the  detailed  energy  transfer 
process  during  a very  short  period  of  time  over  a very  small 
volume  residing  beneath  a growing  vapor  bubble.  Those  are 
the  first  solutions  of  the  kind  according  to  the  literature 
survey.  The  dependence  of  the  temperature  field,  0^,  in  the 
solid  heater  on  the  parameters  (Ja,  Fo,  k,  a)  are  quite 
complicated.  A large  Ja  has  a large  effect  on  0^  due  to  a 
large  amount  of  energy  transferred  from  the  solid  to  the 
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bubble.  A large  Fo  (or  a thin  heater)  implies  a strong 
thermal  interaction  between  the  growing  bubble  and  the 
heater  through  the  microlayer  because  of  the  smaller  thermal 
capacity  of  the  heater.  Therefore,  there  is  a stronger 
energy  flow  within  the  solid  heater  in  the  radial  direction 
to  sustain  the  bubble  growth.  A larger  k also  dictates  a 
stronger  thermal  interaction  because  a larger  gradient  of  6^ 
must  be  maintained  in  order  to  supply  an  adequate  amount  of 
energy  from  the  heater  to  the  bubble  through  the  liquid 
microlayer.  For  a relatively  large  a,  the  heat  conduction 
in  the  radial  direction  is  small  while  in  the  z-direction 
conduction  is  controlled  by  Fo,  not  by  a.  Hence,  for  an 
increasing  a,  the  thermal  diffusion  in  the  z-direction 
becomes  relatively  stronger  when  that  in  the  r-direction 
gets  weaker.  This  results  in  a deeper  penetration  of  the 
solid  temperature  towards  the  bottom  of  the  heater  with 
increasing  a. 

It  is  important  to  keep  in  perspective  that  the  present 

computational  analysis  for  vapor  bubble  growth  is  limited  to 

saturated  pool  boiling  in  the  isolated  bubble  regime. 

Because  the  hydrodynamics  are  not  solved  as  an  integral  part 

of  the  analysis,  it  is  expected  that  the  present  analysis  is 

applicable  only  when  the  Reynolds  number  of  the  bulk  liquid 

motion  induced  by  the  growing  bubble  is  large.  At  low  Jacob 

number  boiling,  the  growth  rate  is  small;  hence  the  Reynolds 

♦ 

number.  Re  = 2R  {t)R  (t) /v  , may  approach  order  unity  which 
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implies  that  both  Reynolds  number  and  Prandtl  number  are 
important  in  determining  C and  . In  fact  it  is  possible 
that  C and/or  may  be  a function  of  time.  Under  these 
extreme  circumstances,  it  is  desirable  to  incorporate  the 
hydrodynamics  with  the  thermal  analyses.  Nevertheless,  the 
present  study  is  quite  useful  because  it  already  covers  a 
wide  range  of  parameters  encountered  in  pool  boiling. 


CHAPTER  4 

A SIMPLIFIED  MODEL  FOR  PREDICTING  VAPOR  BUBBLE  GROWTH 

RATE  IN  HETEROGENEOUS  BOILING 

4 . 1 Background 

In  Chapters  2 and  3 , a detailed  numerical  investigation 
of  the  vapor  bubble  growth  rate  and  associated  thermal 
fields  in  heterogeneous  pool  boiling  has  been  considered. 
Four  dimensionless  parameters  were  identified  to 
characterize  the  vapor  bubble  growth  rate  and  associated 
temperature  fields  of  the  solid  and  liquid.  In  the 
formulation,  a parameter,  C. , that  characterizes  the 
microlayer  wedge  angle  was  determined  by  matching  the  finite 
difference  solution  for  the  vapor  bubble  growth  rate,  R(t) , 
with  the  experimental  data  over  a wide  range  of  boiling 
conditions.  A parametric  investigation  has  elucidated  the 
detailed  energy  transfer  process  among  the  three  phases.  It 
was  shown  that  the  growth  rate  depends  on  the  Jacob  number, 
Ja,  Fourier  number,  Fo,  liquid-to-solid  thermal  conductivity 
ratim,  k,  and  liquid-to-solid  thermal  diffusivity  ratio,  a, 
in  a complex  manner.  No  simple  correlation  could  be 
identified  which  reliably  relates  R(t)  to  these  four 
parameters.  The  vapor  bubble  growth  rate  computation 
described  in  Chapters  2 and  3 is  based  on  two-dimensional 
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axisyininetric  unsteady  heat  conduction  through  a solid  heater 
and  liquid  itiicrolayer . Although  the  bubble  growth 
computation  is  accurate  and  correctly  accounts  for  the 
influence  of  Ja,  Fo,  k,  and  a,  it  is  tedious  and  is  thus 
impractical  to  use  as  a design  tool.  The  need  for  a simple 
bubble  growth  rate  model  which  is  accurate  over  a wide  range 
of  boiling  conditions  is  the  motivation  for  this  work. 

An  energy  balance  on  the  liquid-vapor  interface  in 
conjunction  with  a lumped  thermal  analysis  for  the  solid 
heater  results  in  a first  order  ordinary  differential 
equation  which  is  used  to  calculate  the  vapor  bubble  growth 
rate.  An  empirical  constant,  Cj,  is  introduced  to 
characterize  the  region  of  the  solid  heater  influenced  by 
microlayer  evaporation.  It  is  determined  by  matching  the 
predicted  and  measured  vapor  bubble  growth  rate  data.  It  is 
found  that  a modified  Fourier  number  must  be  introduced  to 
account  for  the  finite  thermal  penetration  depth  when  the 
Fourier  number  is  less  than  one.  Considering  the  complexity 
of  bubble  growth  phenomena,  as  demonstrated  by  Chapters  2 
and  3,  excellent  agreement  is  observed  between  the  predicted 
and  experimental  bubble  growth  data  for  intermediate  and 
moderately  wetting  fluids  in  which  Jacob  number  ranges  from 
0.54  to  1974.  Hence,  for  practical  boiling  problems  of 
interest,  such  as  the  prediction  of  vapor  bubble  departure 
diameters  in  saturated  pool  boiling  Zeng  et  al.  (1993),  the 
simple  model  developed  here  should  be  quite  useful.  As  was 
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the  case  for  the  numerical  model  shown  in  Chapters  2 and  3, 
the  present  correlation  is  limited  to  saturated  pool  boiling 
in  the  isolated  bubble  regime. 


4 . 2 Analytical  Formulation 
4.2.1  Vapor  Bubble  and  Liquid  Microlaver 

Consideration  is  given  to  a vapor  bubble  growing  from  a 
solid  heating  surface  into  a saturated  pool  of  liquid.  The 
bubble  shape  is  assumed  to  consist  of  a truncated  sphere  on 
a wedge-shaped  microlayer  as  was  assumed  in  Chapter  2 , and 
shown  in  Fig.  2.1.  The  apex  of  the  microlayer  wedge  is 
assumed  to  be  at  r=0  as  shown  in  Fig.  2.1.  This  assumption 
alleviates  the  requirement  to  specify  the  movement  of  the 
three-phase  contact  point  and  thus  greatly  simplifies  the 
analysis.  As  long  as  the  vapor-solid  contact  area  is  small 
comparing  with  the  liquid-solid  areas  in  the  wedge,  this 
assumption  will  not  cause  significant  error  in  predicting 
the  heat  transfer.  In  the  parametric  investigation 
conducted  in  Chapter  3,  it  was  established  that  the 
temperature  profile  in  the  liquid  microlayer  is  practically 
linear.  Hence,  the  liquid  temperature  in  the  microlayer  is 
approximately 


\ 


1--^ 


(4.1) 


where  AT^^^  (r , t)  =Tg (r , z=0,  t) -T^^^,  T^^^  is  the  saturation 
temperature,  Tg(r,z,t)  is  the  solid  heater  temperature,  L(r) 
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is  the  local  microlayer  thickness,  and  t denotes  time. 

Thus,  the  heat  flux  at  the  microlayer  liquid-vapor  interface 
is. 


z=Z(/) 


for  <j)<l  , 


(4.2) 


where  it  has  been  assumed  that  the  microlayer  wedge  angle, 
0,  is  very  small  as  shown  in  Fig.  2.1.  It  is  also  assumed 
that  the  evaporation  of  the  microlayer  dominates  the  heat 
transfer  to  the  vapor  bubble  during  the  bubble  growth 
period.  Thus,  the  change  in  the  bubble  volume  can  be 
obtained  through  an  energy  balance  on  the  liquid-vapor 
interface. 


dA 

z=£{r) 


(4.3) 


where  h^  is  the  latent  heat,  is  the  vapor  density, 
is  the  microlayer  surface  area.  Denoting  Rt,(t)  the 
microlayer  base  radius,  R(t)  as  the  bubble  radius,  as 
in  Fig.  2.1,  the  bubble  volume  can  be  expressed  as 

V, = ( 1 -|(i  ) 


and  A 


shown 


(4.4) 


where  C=R|j(t) /R(t)  <1  is  the  bubble  shape  factor,  as  was 
shown  in  Chapter  2.  Since  dA«27rrdr  and  L(r)ar0,  (4.3)  may 
be  expressed  as. 
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dR'^  Cl 


for  <f><l . 


dt  f{C)<i>  ?Jifg 


(4.5) 


In  the  above,  AT  (t)  is  the  superheat  on  the  boiling 

SdL 

surface  averaged  over  the  area  beneath  the  microlayer  and  is 
obtained  from  an  energy  balance  within  the  solid  heater.  If 
the  wall  superheat  was  constant  during  the  bubble  growth 
period,  AT^^=AT^^^,  and  equation  (4.5)  may  be  used  to  solve 
for  the  growth  rate. 


where  the  subscript  c denotes  the  solution  for  a constant 
wall  superheat,  is  the  liquid  thermal  diffusivity,  and  Ja 
is  the  Jacob  number  based  on  the  initial  wall  superheat, 

AT  _Q.  Typically,  this  superheat  is  either  specified  or 
determined  by  considering  the  total  heat  transfer  from  the 
solid  heater.  The  microlayer  wedge  angle  cp  was  generalized 
in  chapter  2 based  on  the  analysis  of  Cooper  & Lloyd  (1969) 
as 


(4.6) 


4>= 

C^Ja 


(4.7) 


where  Pr^  is  the  liquid  Prandtl  number,  and  is  an 
empirical  parameter  (see  Eq.  (4.17))  determined  (Chapter  2) 
by  matching  the  predicted  and  measured  growth  rate.  The 
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case  for  time  varying  wall  superheat  is  considered  next. 

4.2.2  Lumped  Analysis  for  the  Solid  Heater 

Most  of  the  vapor  bubble  growth  data  available  in  the 
literature  were  obtained  using  thin  heaters  in  which  the 
temperature  on  the  back  side  is  almost  equal  to  that  on  the 
boiling  side.  Consideration  is  given  to  thin  heaters  first. 
Under  such  circumstances,  a "lumped  capacitance"  approach 
may  be  used  to  formulate  the  energy  balance  for  the  solid 
heater  with  thickness  H as 


V p C 

c r c / 


dAT 


sat 


= q'"  Vs- 


ps 


dt 


evap 


(4.8) 


where  V is  the  effective  volume  in  the  solid  influenced  by 
the  microlayer  evaporation,  AT^^^t)  is  understood  to  be  the 
volume  averaged  superheat  in  V^,  is  the  internal  heat 

generation  per  unit  volume  (or  = q^^  /H  in  which  q^^  is 
the  heat  flux  supplied  from  the  bottom  of  the  heater)  , 
is  the  rate  of  heat  removal  from  the  heating  surface  due  to 
microlayer  evaporation.  The  subscript  s denotes  that  the 
quantities  are  for  the  solid  heater.  It  is  necessary  to 
specify  an  effective  heater  volume,  Vg,  from  which  energy  is 
extracted  in  order  to  satisfy  the  energy  demand  from  bubble 
growth.  Here  it  is  assumed  that  Vg=7r(C2R)^H  in  which  C2  is  a 
constant  to  be  determined.  Combining  equations  (4.3)  and 
(4.5),  and  noting  that  q is  the  right-hand-side  of 
equations  (4.3),  (4.8)  becomes 


(4.9) 


dLT. 


sat 


q!"  \f{C)^JifgdR 

Cip,C^Hdt' 


4.2.3  Non-dimensionalization 

Following  Chapter  2,  the  dimensionless  variables  (see 
nomenclature)  are  introduced  to  recast  equations  (4.5)  and 


(4.9)  as 


and 


dr 


dr 


Q-B^ 

R 


where 


Q = 


q'"t 


^^satO^  s^p,s 


2k  if{C)CFo'^r^ 

\ aJa4> 


(4.10) 

(4.11) 

(4.12) 


In  equation  (4.12),  Fo  is  the  Fourier  number  which  is 
defined  previously.  Eliminating  6 from  (4.10)  and  (4.11),  a 
first  order  ordinary  differential  equation  for  R{z)  results. 


dR 

dr 


^2BR=  QT^\^2BR{Qi) 


(4.13) 


where  J?(0)  is  the  initial  bubble  radius  which  is  physically 
very  small  and  is  mathematically  set  to  zero  in  this  study. 
Equation  (4.13)  is  easily  solved  using  a Runge-Kutta  method. 
When  Q is  negligibly  small  or  q^^^  is  negligible  compared  to 
the  heat  removed  through  microlayer  evaporation,  an  exact 
solution  to  equation  (4.13)  is  obtained. 
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T 


^ ’ B 


(4.14) 


For  small  B or  small  R , a Taylor  series  expansion  for  the 
logarithmic  term  results  in 


R{t)  ~ + H.O.T. 


(4.15) 

which  is  the  dimensionless  form  of  Rj.(t)  given  by  equation 
(4.6).  Since  we  are  interested  in  R{t)~  0(1)  and  in  general 
B is  not  small  unless  Ja  is  very  large,  equation  (4.6)  or 
(4.15)  is  not  a satisfactory  solution  for  t ~ 0(1).  This 
is  the  reason  the  solution  based  on  the  constant  wall 
superheat  is  unsatisfactory  when  compared  with  experimental 
data  over  a wide  range  of  Jacob  number.  It  is  also  noted 
that  for  a thick  heater,  B-^0  as  Fo-»0,  and  equation  (4.14) 
again  degenerates  to  the  constant  wall-superheat  solution. 
Hence  it  fails  to  predict  the  actual  bubble  growth  rate. 
From  the  previous  two-dimensional  computation  (Chapters  2 
and  3),  it  is  known  that  for  a small  Fo,  the  penetration 
depth  of  the  thermal  field  in  the  solid  heater  is  actually 
somewhat  limited  prior  to  departure.  Hence  it  is 
anticipated  that  when  Fo  is  less  than  one,  a modified, 
effective  Fourier  number  should  be  used. 
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4 . 3 Results  and  Discussion 

4.3.1  Results  for  Fo>l;  Determination  of  the  Empirical 
Constant 

The  empirical  constant  Cj  mainly  characterizes  the 
region  of  the  heating  surface  influenced  by  microlayer 
evaporation;  it  also  partially  accounts  for  the  non- 
uniformity in  the  heater  thermal  field  beneath  the 
microlayer.  It  is  determined  by  matching  the  predicted 
growth  rates,  computed  from  equation  (4.13),  with  the 
experimental  growth  rate  data.  The  two  empirical  constants, 
C and  , have  been  determined  previously  (see  chapter  2)  as 

C=  ((0.4134 /a  (1  -0.1 

(4.16) 

and 

q Prf  = 0.00525 


There  are  32  data  listed  in  Chapter  2.  Nine  of  them  have 
Fo<l.  Based  on  the  remaining  23  data  for  Fo>l, 
determined  by  best  fitting  the  predicted  growth  rate  with 
the  measured  ones.  The  results  for  Cj  are  given  in  Table 
4.1.  Based  on  the  results  for  ^ correlation  is 
developed  as 


q = 1.26  (0.73 + 


(4.18) 


This  gives  a relative  deviation 


Table  4 . 1 The  empirical  constant  C,  determined  from  various  sources 
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Table  4.1  (continued) 
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1 23  ic  -q*  I 

r.d  = —Y  ' ~ 0.085 

nk  a, 

where  C^n  is  the  value  given  in  Table  4.1  which  is 
determined  by  matching  the  measured  and  predicted  growth 
rates.  As  was  mentioned  before,  C2  is  a constant  which 
characterizes  the  effective  volume  in  the  solid  influenced 
by  the  microlayer  evaporation.  It  is  observed  from  the  data 
that  Cj  is  mainly  correlated  with  Fo  and  Ja:  C2  increases 
with  increasing  Fo  and  decreasing  Ja.  As  the  heater  gets 
thinner  (which  increases  Fo) , the  energy  transfer  has  to 
rely  more  on  the  radial  conduction  than  the  vertical 
conduction  in  order  to  sustain  the  vapor  bubble  growth. 

This  leads  to  a larger  spreading  of  the  region  of  influence 
by  the  microlayer.  The  effect  of  increasing  (or  decreasing) 
Ja  on  the  radial  thermal  diffusion  has  been  discussed  in 
detail  in  Chapter  3;  it  was  observed  that  as  Ja  decreases 
the  radial  conduction  is  enhanced.  This  is  consistent  with 
the  trend  observed  for  the  C2  dependence  on  Ja. 

To  validate  the  present  simplified  model  based  on  the 
lumped  analysis,  comparisons  are  first  made  between  the 
predicted  R(t) , which  includes  the  present  prediction 
(equation  4.13)  and  the  computational  prediction  from 
Chapter  2,  and  the  measured  R(t)  for  Fo>l.  The  results  are 
shown  in  Figures  4.1a-c.  It  is  noted  that  for  the  predicted 
growth  rate,  R(t) , shown  in  Figure  4.1,  C2  was  calculated 
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Figure  4.1a  Comparision  of  the  bubble  growth  rate  R(t) 
between  the  present  model,  two-dimensional  model,  and  the 
experimental  data  of  Cole  and  Shulman  (1966)  using  water  as 
boiling  liquid. 
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Figure  4.1b  Comparison  of  the  bubble  growth  rate  R(t)  between 
the  present  model,  two-dimensional  model,  and  the  experimental 
data  of  Cole  and  Shulman  (1966)  using  methanol  as  boiling 
liquid. 


R(t)  ( 
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Figure  4.1c  Comparison  of  the  bubble  growth  rate  R(t)  between 
the  present  model,  two-dimensional  model,  and  the  data  of 
Staniszewski  (1959)  using  water  as  boiling  liquid. 
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using  equqtion  (4.18)  (values  for  C2  given  in  Table  1 were 
not  used) . However,  is  taken  from  the  values  listed  in 
Table  4.1.  In  this  way,  the  difference  between  the 
predictions  based  on  the  present  model  and  that  based  on  the 
two-dimensional  model  is  only  due  to  the  modeling  of  the 
lumped  system  and  the  correlation  for  C2,  not  . In  Figure 
4.1a,  the  data  are  from  Cole  & Shulman  (1966)  in  which  Ja 
ranges  from  intermediate  to  high  values,  water  is  the 
boiling  liquid,  and  the  heater  material  is  a polished 
zirconium  ribbon.  Excellent  agreement  for  R(t)  is  observed 
except  for  data  No.  12.  In  the  case  of  data  No.  12,  the 
relative  deviation  between  C2  based  on  the  correlation 
equation  (4.18)  and  the  best  fit  value  listed  in  Table  4.1 
is  about  17%,  which  is  the  reason  for  the  poor  agreement. 

The  comparison  with  data  No.  12  is  included  here  because  the 
agreement  is  the  worst  among  the  32  sets  of  data.  Figure 
4 . lb  compares  the  present  prediction  with  the  data  of  Cole 
and  Shulman  (1966)  in  which  Ja  is  considered  in  the 
intermediate  range  and  methanol  is  the  boiling  liquid.  Very 
good  agreement  is  again  observed.  Fig.  4.1c  compares  the 
present  prediction  with  the  data  of  Staniszewski  (1959)  in 
which  water  is  the  boiling  liquid,  the  heater  material  is 
copper  and  Ja  is  considered  low.  Good  agreement  is 
observed  for  most  of  the  growth  period. 

4.3.2  Results  for  Fo<l;  Modified  Fourier  Number  ( Fo^) 

The  bubble  growth  computation  of  chapter  3 indicates 
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that  the  bubble  growth  rates  i?  (r)  for  Fo=0.01  and  for  Fo=l 
are  almost  the  same;  but  -R  ('c)  for  Fo=l  are  larger  than  for 
Fo=10  when  other  parameters  are  held  the  same.  The  reason 
is  as  follows:  once  the  thickness  of  the  heater  exceeds  some 
value,  the  thermal  field  in  the  bottom  portion  of  the  heater 
is  unaffected  by  the  growth  of  the  vapor  bubble  during  the 
entire  growth  period,  r<l.  For  such  cases,  H in  the  present 
lumped  analysis  must  be  modified  to  reflect  the  effective 
thermal  penetration  depth  which  is,  unfortunately,  time 
dependent.  To  alleviate  this  difficulty,  we  use  a modified 
Fourier  number,  Fo^j.  When  using  equation  (4.18)  to  evaluate 
C2,  Fo  is  replaced  by  Fo^  in  which 


Fo 


= < 


1.29  - 0.29exp(l 


Fo>\ 

Fo<\ 


(4.20) 


and  Fo  also  replaces  Fo  in  equations  (4.12)  and  (4.13).  It 
is  noted  that  in  the  limit  Fo->l,  Fo^-^1;  as  Fo-*^0,  Fo^-^0.5. 

For  the  remaining  nine  data  with  Fo<l,  which  can  be 
identified  in  Table  4.1,  the  present  predictions  for  R(t) 

agree  well  with  the  measured  ones.  Figure  4 . Id  compares 

the  predicted  R(t)  and  the  measured  R(t)  for  Fo<l  in  the 
high  Jacob  number  range.  The  data  are  from  Van  Stralen  et 

al.  (1975)  with  water  as  the  boiling  liquid  on  a copper 

heater.  As  seen,  the  predicted  bubble  growth  rates  shown  in 
Figure  4 . Id  compare  well  with  the  measured  data. 


R(t)  ( 
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Figure  4 . Id  Comparison  of  the  bubble  growth  rate  R(t)  between 
the  present  model,  two-dimensiional  model,  and  the  data  of  Van 
Stralen  et  al.  (1975)  using  water  as  boiling  liquid  for  Fo<l. 
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4 . 4 Summary 

In  this  study,  a simplified  model  is  formulated  to 
predict  the  growth  rate  of  vapor  bubbles  in  saturated 
heterogeneous  pool  boiling  for  the  isolated  bubble  regime 
based  on  the  same  physical  mechanisms  identified  in  Chapters 
2 and  3.  A first  order  ordinary  differential  equation  is 
obtained,  which  when  numerically  solved,  yields  the  vapor 
bubble  growth  rate.  The  predicted  bubble  growth  rates 
compare  well  with  the  reported  experimental  data  in  which 
the  working  fluids  are  intermediate  and  moderately  wetting 
and  the  Jacob  number  ranges  from  0.52  to  1974.  A 
correlation  for  Cj,  which  characterizes  the  region  of  the 
heating  surface  influenced  by  the  microlayer  evaporation,  is 
presented.  For  thick  heaters  with  Fo<l,  a modified  Fourier 
number  is  introduced  to  account  for  the  finite  thermal 
penetration  depth  which  allows  the  vapor  bubble  growth  rate 
to  be  accurately  predicted.  The  present  bubble  growth  model 
provides  a simple  and  apparently  reliable  alternative  to  the 
intensive  two  dimensional  finite  difference  computation 
presented  in  Chapters  2 and  3 . Due  to  the  empiricism 
involved  in  the  model  development,  it  must,  nevertheless,  be 


used  with  caution. 


CHAPTER  5 

VAPOR  BUBBLE  GROWTH  IN  SUBCOOLED 
HETEROGENEOUS  BOILING 

5 . 1 Background 

In  previous  investigations  discussed  in  Chapters  2 and 
3,  detailed  theoretical  and  numerical  analyses  of  vapor 
bubble  growth  in  saturated  heterogeneous  pool  boiling  have 
been  studied  over  a wide  range  of  boiling  conditions.  In 
the  formulation,  the  simultaneous  coupling  of  the  thermal 
fields  in  a thin  liquid  microlayer  and  in  the  solid  heater 
takes  into  account  the  continuously  decreasing  heating 
surface  temperature  beneath  the  growing  bubble.  In  the 
parametric  investigation,  four  dimensionless  parameters  were 
identified  to  characterize  the  vapor  bubble  growth  rate  and 
associated  temperature  fields  in  the  solid  and  liquid 
microlayer.  The  numerical  analysis  has  elucidated  the 
detailed  energy  transfer  process  for  the  solid  thermal 
field.  In  saturated  pool  boiling,  the  bulk  liquid 
surrounding  the  bubble  cap  is  nearly  saturated,  and  thus  the 
heat  and  mass  transfer  in  this  region  is  small.  In 
subcooled  boiling,  the  unsteady  heat  transfer  between  the 
vapor  bubble  dome  and  the  bulk  liquid  is  no  longer 
negligible  due  to  the  effect  of  subcooling.  During  the 
early  stages  of  growth,  the  vapor  bubble  dome  obtains  energy 
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from  the  surroundings  because  a large  portion  of  the  bubble 
surface  area  is  within  the  thermal  boundary  layer.  Thus, 
the  overall  energy  for  the  bubble  growth  comes  from  the 
liquid  microlayer  and  the  vapor  bubble  dome.  As  the  bubble 
continues  to  grow,  a large  portion  of  the  bubble  dome 
protrudes  the  thermal  boundary  layer.  At  the  same  instant, 
the  energy  into  the  bubble  from  the  microlayer  is  smaller 
than  that  lost  from  the  dome  to  the  surrounding  subcooled 
liquid,  and  the  bubble  starts  to  collapse.  In  this  study, 
the  vapor  bubble  growth  process  is  considered  to  be 
controlled  by  the  unsteady  energy  transfer  among  four 
regions:  a wedge  shaped  liquid  microlayer,  a vapor  bubble 
whose  shape  is  assumed  to  be  a truncated  sphere,  a subcooled 
bulk  liquid,  and  a solid  heater.  A linear  temperature 
profile  across  the  liquid  microlayer  is  assumed  based  on  the 
investigation  in  Chapter  3 . The  unsteady  energy  equation  in 
the  subcooled  liquid  which  is  assumed  to  be  axisymmetric 
with  respect  to  the  axis  normal  to  the  heating  surface,  is 
strongly  coupled  with  the  momentum  equation.  The  momentum 
equation  is  approximated  by  the  Laplace  equation  by  assuming 
an  inviscid  flow.  Good  agreement  between  the  predicted  and 
measured  bubble  radius  is  achieved.  The  present  study 
contributes  to  the  basic  understanding  of  the  bubble  growth 
mechanisms  during  subcooled  nucleate  boiling. 
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5 . 2 Theoretical  Formulation 
5.2.1  Vapor  Bubble  and  Liquid  Microlaver 

Consideration  is  given  to  an  isolated  vapor  bubble 
growing  from  a solid  heating  surface  into  a subcooled  liquid 
pool  as  shown  in  Figure  5.1.  The  bubble  shape  is  assumed  to 
be  a truncated  sphere  on  a wedge-shaped  microlayer.  This 
assumption  alleviates  the  requirement  to  specify  the 
movement  of  the  three-phase  contact  point  and  thus  greatly 
simplifies  the  analysis.  In  the  parametric  investigation 
for  saturated  boiling  (Chapter  2),  it  was  established  that 
the  temperature  profile  in  the  liquid  microlayer  is 
practically  linear  and  the  liquid  temperature  in  the 
microlayer  is  given  by 


T,{r,z,t)  = 7;^,+ 


(5.1) 


where  AT^^Jr,  t)  =Tg(r,  z=0,  t) L(r)=r0  is  the  local 
microlayer  thickness,  and  <p  is  the  microlayer  wedge  angle 
which  is  generally  very  small.  During  the  bubble  growth  and 
collapse,  the  change  in  the  bubble  volume  can  be  obtained 
through  an  energy  balance  on  the  microlayer  liquid-vapor 
interface  and  the  subcooled  bulk  liquid-vapor  interface. 


dt 


-L 


-k. 


m 


dn 


L 


dT, 


bi 


dn 


R'=R{t) 


dA^  , for  <^<1  , 


(5.2) 


where  T|^^  is  the  temperature  of  the  subcooled  liquid,  is 
the  surface  area  of  wedge  on  the  microlayer  liquid-vapor 


Figure  5.1  Sketch  for  the  growing  bubble,  liquid 
microlayer,  subcooled  bulk  liquid,  and  heating  solid. 
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interface,  is  the  surface  area  of  the  vapor  bubble  dome 
on  the  subcooled  bulk  liquid-vapor  interface,  b/bn  is 
differentiation  along  the  outward-drawn  normal  at  the 
boundary  surface,  and  is  the  spherical  coordinate  in  the 
radial  direction  attached  to  the  moving  bubble.  Equation 
(5.2)  simply  states  that  the  energy  conducted  from  the  solid 
to  the  bubble  through  the  microlayer  and  the  energy  transfer 
between  the  vapor  bubble  and  subcooled  bulk  liquid  cause  the 
vapor  bubble  to  grow  or  condense.  Denoting  R^Ct)  as  the 
microlayer  base  radius  and  R(t)  as  the  bubble  radius,  the 
bubble  volume  can  be  expressed  as 

n = ( 1 -|(  1 -^1  f f ) 

where  C(t)=Rjj(t) /R(t)<l  is  the  bubble  shape  factor,  as  was 
mentioned  in  chapters  2 and  4.  In  the  limit  as  the  bubble 
shape  approaches  a hemisphere,  C(t)-+1,  and  as  the  bubble 
shape  approaches  a sphere,  C(t)^0.  Based  on  experimental 
evidence  by  Cochran  and  Aydelott  (1966)  and  Ellion  (1954) , 
the  bubble  shape  in  subcooled  pool  boiling  changes  with  time 
from  a hemispherical  dome  in  the  early  stage  towards  a 
spherical  dome.  Thus,  C varies  with  time. 

5.2.2  Subcooled  Bulk  Liquid 

By  assuming  axisymmetry  for  the  temperature  and 
velocity  fields,  the  unsteady  energy  equation  for  the 
subcooled  bulk  liquid  in  cylindrical  coordinates  (r,z)  is 
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(5.4) 


where  is  the  thermal  diffusivity,  (u^,u^)  are  the  liquid 
velocity  components  in  the  (r,z)  directions.  The  initial 
condition  specified  for  the  liquid  thermal  field  is  based  on 
Marcus  and  Dropkin's  (1965)  investigation  of  the  thermal 
boundary  layer  above  a copper  surface  in  saturated  nucleate 
pool  boiling  of  water: 


where  is  the  average  temperature  of  the  boiling  surface, 
w and  X are  experimentally  determined  constants  based  on  the 
measurements  of  the  subcooled  liquid  temperature, 

Z>=(1 -0.57 0.57^  , and  S is  the  extrapolated  superheat-layer 
thickness.  It  is  noted  that  w=l  in  saturated  nucleate  pool 
boiling.  In  subcooled  pool  boiling,  equation  (5.5)  is 
adopted  by  noting  that  the  value  of  w,  which  can  be  obtained 
from  the  data  of  Wiebe  and  Judd  (1971)  may  be  a function  of 
heat  flux  and  degree  of  subcooling  and  varies  from  1 to  1.75 
for  the  data  considered  here.  X is  set  to  be  1.2  in  the 
numerical  computation.  The  extrapolated  superheat-layer 
thickness,  5,  is  estimated  by  Hsu  (1962)  using  the 
expression  given  which  was  obtained  assuming  the  thermal 


-<0.57 

8 


->0.57 

8 


(5.5) 


layer  is  governed  by  one-dimensional  unsteady  heat 
conduction. 
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'^subt^ 


IACmIAQ 


+ 


(5.6) 


where  Cj  is  equal  to  1.6,  A is  equal  to  2aTgg^/p^hfg,  and  a is 
the  surface  tension.  Hsu  has  demonstrated  the  usefulness  of 
equation  (5.6)  for  the  prediction  of  the  incipience  point 
over  a range  of  subcoolinq  and  pressure.  As  suggested  by 
Wiebe  and  Judd  (1971)  equation  (5.6)  is  used  in  this  study 
to  estimate  the  extrapolated  superheat-layer  thickness  S . 

Since  there  is  no  strong  mean  flow  over  the  bubble,  the 
flow  induced  by  the  growth  of  the  bubble  is  mainly  of 
inviscid  nature.  Thus  the  liquid  velocity  field  may  be 
determined  by  solving  the  Laplace  equation  V^<I>=0  for  the 
velocity  potential  $ and  the  velocities  are  given  by. 


= 


dr 


and  u^r,d) 


dz 


(5.7) 


Since  the  bubble  is  a growing  truncated  sphere,  the  velocity 
field  needs  to  be  solved  at  each  time  step.  If  the  vapor 
bubble  is  exactly  hemispherical,  which  is  the  case  during 


the  early  growth  stage,  the  liquid  velocities  are 


=R{t) 


^ sin^(i|;)  and  uj^r,\if)  = R(i)  ^ cos^(i}f)  , 

r^  r^  (5.8) 


dR(t) 


where  R(t) 


dt 


and  \{/  is  measured  from  the  axis  of 
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symmetry,  i.e.,  z-axis.  The  general  boundary  conditions  for 
$ will  be  given  in  the  transformed  coordinates  when  a non- 
hemispherical  bubble  is  considered. 

5.2.3  Solid  Heater 

The  energy  equation  for  the  solid  is 


dJ^ 

dt 


dz^ 


+ —q'". 


k. 


0 < r<  00,  -H<z<Q 


(5.9) 


subject  to  the  initial  condition 
and  the  boundary  conditions 


at  t=L<~t 

U C 


dr 


at  i*=0  and  r-^oo, 


(5.11) 


-^  = -qiilk,  atz=-H, 

dz  ^ ^ 


(5.12) 


at  z=0,  Q<r<Rij{t), 


(5.13) 


-7;-^  = 7? ( T; - T;„ J at  z=Q, 

(5.14) 

where  and  are  the  thermal  diffusivity  and  conductivity 

O w 

of  the  solid,  t^.  is  the  bubble  collapse  time  scale,  Tg^^^Q  is 
the  subcooled  bulk  liquid  temperature  at  infinity,  h is  the 
average  convective  heat  transfer  coefficient  in  the 
subcooled  bulk  liquid,  q^^  is  the  heat  flux  (on  the 
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average)  , and  is  the  volummetric  heating  (on  the 

average) . 

5.2.4  Non-dimensionalization  and  Grid  Generation 

To  obtain  a better  understanding  of  the  physics  of  the 
problem,  the  following  dimensionless  variables  are 
introduced, 


T = 


r. 


r 


R(t) 


z.  = 


z 


H 


T - T T - T 

\ *' J /I  /\  bC  subO 


- m 


(5.15) 


at: 


satO 


T -T 

sat  subO 


where  R.  (t.)  is  defined  in  Chapter  2.  Following  Chapter  2, 
by  using  equation  (5.15)  and  applying  coordinate 
transformations , 


(l-5;tan-i[(l-0tan(l/5;)]}  for  0<C<1, 

(5.16) 

= -1  + 5'^tan“^[(l-£)tan(l/5y)]  for  0<£<1, 

(5.17) 

and  T=a^,  equation  (5.9)  becomes 

1 ac  _ 1 Idlj  I 0^ 

2a  da  ^ dC  dJ^  aJaCfij)  R^{t)  ^5  ' "'SC  drjdr^ 

ds\dsdzjdz^ 


where  a and  k are  defined  in  the  nomenclature,  is  the 
radial  position  at  which  equation  (5.9)  is  applied,  the 
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constant  S^.  determines  the  percentage  of  the  grid  number 
allocated  to  the  region  of  / and  is  a constant 

characterizing  the  stretching  in  the  z-direction. 

For  the  subcooled  bulk  liquid,  the  energy  equation  in 
dimensionless  form  is 


To  solve  the  Laplace  equation  for  $ and  energy  equation  for 
6yf,  an  algebraic  grid  generation  technique  is  used  for  the 
bulk  liquid  outside  the  growing  bubble.  Figure  5.2  shows 
the  computational  grids  for  the  numerical  simulation.  The 
clustering  of  grids  near  the  bubble  surface  is  necessary  to 
resolve  large  temperature  gradients.  The  and  within 
the  physical  domain  are  calculated  as 


= A/)  sin(Y(v))  > 

Zteiij)  = A/)  cos(y(v))  + z^, 


(5.20) 


where 


(5.21) 
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Figure  5.2  The  computational  grids  for  the  numerical 
simulation. 
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XV 

is  the  radial  position  at  infinity  and  the  constant 
determines  the  percentage  of  the  grid  number  near  the  bubble 
surface.  Figure  5.3  illustrates  the  grid  parameters. 

The  transformed  energy  and  Laplace  equations  for  the 
subcooled  bulk  liquid  are 


(5.22) 


and 


y^l  a +2o + c — + a — + e — 1+  — — -z,„  — 1 = U 


a^2  ^2  dTj  Qi!  r^\ 


'Hi 


dx] 


(5.23) 


where  /=  — — _ 


^Hf  ^H^  ^Hf 


which  is  the  Jacobian, 


^ ^H^~^  ^H^’  ^ ^hJh.  ^ ^H J^H{  ^ ^Hf  ^ ^H{’ 

TJTJ  1?C^C  C C 


CC  Wl 


The  temperature  boundary  conditions  are  along  the 
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Figure  5 . 3 Grid  parameters 
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vapor-liquid  interface,  and  the  temperature  gradient  is  zero 

along  the  line  of  symmetry.  The  heat  flux  at  the  solid 

heater  is  equal  to  the  heat  flux  to  the  subcooled  bulk 

liquid  along  the  wall.  The  far-field  boundary  condition  for 

is  approximated  by  zero  heat  flux.  The  boundary 

conditions  for  velocity  potential  are  zero  velocity 

potential  gradient  along  the  wall  and  the  line  of  symmetry. 

The  far-field  boundary  condition  is  approximated  by  a zero 

flux  condition.  The  boundary  condition  at  the  vapor-liquid 

interface  is  — = i?  + Acos  il; , where  h = is  the  bubble 

dn  dt 

vertical  velocity  due  to  its  moving  center,  and  h(t)  is  the 
distance  from  the  wall  to  the  bubble  center.  Equation 
(5.22)  is  solved  using  the  ADI  method  (Anderson  et  al., 

1984)  and  equation  (5.23)  is  solved  using  SOR  by  line 
(Anderson  et  al.,  1984). 

5 . 3 The  Asymptotic  Analysis  of  6^^ 

To  gain  a basic  understanding  of  the  behavior  of 
0|^(t,r,z)  associated  with  a rapidly  growing  vapor  bubble,  an 
asymptotic  analysis  for  is  presented  in  this  study.  For 
simplicity,  we  consider  the  growth  of  a hemispherical 
bubble.  As  mentioned  early,  during  the  earlier  stage  of  the 
bubble  growth,  the  bubble  shape  is  indeed  hemispherical. 

Due  to  the  rapid  bubble  growth  and  the  limited  extent  of  the 
heat  diffusion  effect,  there  exists  a thin  thermal  boundary 
layer  adjacent  to  the  bubble  surface.  Thus  an  asymptotic 
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analysis  can  be  used  to  describe  the  complex  phenomena. 

Using  the  dimensionless  variables  in  equation  (5.15) 

“/ 

and  defining  a new  dimensionless  variable  as  R = , the 

energy  equation  for  the  subcooled  liquid  in  spherical 
coordinates  is, 


R 

t^k 


a, 


/ 


2 » 


\ 


/2 


+ 


de 


r'  dR' ) 


1 1 ^ 

^/2  aijj 


(5.24) 


\ 

During  the  early  stage  of  growth,  the  bubble  growth 
rate  is  high  and  expands  rapidly  so  that 


(5.25) 


Thus,  the  solution  for  equation  (5.24)  includes  the  outer 
approximation  and  inner  approximation  (boundary  layer 
solution) . To  obtain  the  outer  solution  for  equation 
(5.24),  the  diffusion  term  in  equation  (5.24)  is  neglected. 
This  results  in,  to  the  leading  order, 

out  out 

R dt  jr^/2  (5.26) 

The  general  solution  of  equation  (5.26)  is 


(5.27) 


Ill 


where  F is  an  arbitrary  function  and  it  is  determined  from 

the  initial  condition.  It  is  noted  that  the  solution  for 

is  described  by  -1)^^^  = constant  along  the 

characteristic  curve  . Due  to  the  discontinuity  in  of 

dz 

the  initial  condition,  equation  (5.5),  we  use  the  following 
continuous  function  instead  of  equation  (5.5),  in  discussing 
the  analytical  solution  for  6^^, 


-T 


subO 


-T 


subO 


(5.28) 


This  function  matches  the  linear  part  of  equation  (5.5)  and 
the  computational  results  indicate  that  there  only  exists 
about  3~4%  difference  in  the  computed  bubble  radius  R(t) 
when  equation  (5.28)  is  used.  To  determine  F,  the  initial 
condition  given  by  equation  (5.28)  is  employed. 


6/ 


out 

be 


(5.29) 


where 

T - T 

^ ^ yy  ^ subO 

subO  rp  _ rp  f 

sat  subO 

and  Rq  is  the  initial  bubble  radius  at  t=tg<tj.  and  Rg  can  be 
easily  obtained  from  equation  (5.2).  As  long  as  the  bubble 
growth  rate  is  high,  i.e.  A>1,  equation  (5.29)  is  not  only 
an  accurate  outer  solution  for  a hemispherical  bubble,  but 
also  is  a good  approximation  for  the  far  field  boundary 
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condition  for  6^^  in  equation  (5.22)  for  non-hemispherical 
bubbles . 

The  outer  solution  automatically  satisfies  the 
condition  for  6^^  as  R^-*<x>.  However,  on  the  bubble  surface, 
i?=l,  cannot  be  satisfied  by  The  smooth 

transition  of  6^^  from  6^“^  (R^  = 1)  to  R^=l  calls  for 

an  inner  (or  boundary  layer)  solution.  In  this  thin 
boundary  layer,  the  effect  of  heat  conduction  is  no  longer 
negligible  and  must  be  properly  accounted  for.  For  large  A, 
a boundary  layer  coordinate  X is  introduced 


-I 

S\A) 


(5.30) 


where  5*  (A)  is  the  scale  of  the  boundary  layer  which  is 
determined  by  balancing  the  convection  term  with  the 
conduction  term  in  the  boundary  layer.  Substituting  equation 
(5.30)  into  equation  (5.24)  results  in 
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Since 
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equation  (5.31)  becomes 
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in  the  limit  5*(A)-»'0  as  A-«».  The  balance  between  the 


convection  term  and  the  diffusion  term  requires  5 =A 
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Thus,  equation  (5.33)  results  in 
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subject  to  the  boundary  conditions 
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at 
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For  T<^1,  R(t)  «t^^^  and 
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~T.  Thus,  the  LHS  of  equation 


(5.34)  is  small  and  equation  (5.34)  reduces  to 
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The  solution  for  equation  (5.37)  is 


The  uniformly  valid  solution  for  6^^  for  t<^1  is 

-^u-  since  = we  have 

««  = 1 * («„«  - 1)  X) 


\ i^cosil; 

+ ^sub0^^v[-^ ^ 


(5.38) 


(5.39) 


Equation  (5.39)  is  an  approximate  asymptotic  solution  for  0^ 

• 

for  As  RIR  increases,  equation  (5.39)  becomes 

increasinqly  inaccurate  due  to  the  neqlect  of  the  LHS  in 
equation  (5.34). 


5 . 4 Results  and  Discussion 

The  shape  parameter  C(t)  mainly  describes  how  the  shape 
of  the  vapor  bubble  changes  during  the  bubble  growth  and 
collapse.  Due  to  a lack  of  experimental  information,  the 
bubble  shape  factor  in  this  study  is  assumed  to  vary  as 
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where  t is  the  time  when  the  bubble  shape  begins  to  change 
from  a hemisphere  to  a sphere.  The  determination  for  t^ 
employs  a dimensionless  parameter  fi  = ^ which 

characterizes  the  relative  growth  rate.  It  is  noted  that  /Xg 
(i.e.  t=t  ) is  determined  by  matching  the  computed  growth 
rates  with  the  experimental  data  in  order  to  achieve  the 
best  fit.  In  this  study  ranges  from  0.23  to  0.48. 

Two  empirical  parameters,  0 and  /Xg,  need  to  be 

specified  for  the  present  computation.  Through  the 

comparison  of  predicted  growth  rates  with  experimental 

results,  0 and  /x„  can  be  determined  via  the  best  fit.  From 

the  numerical  results  there  exists  only  one  combination  of 

0 and  /x„  for  the  best  agreement  with  experimental  data. 

Although  the  present  computational  model  seems  crude  due  to 

the  uncertainty  in  the  two  empirical  constants,  the  main 

contribution  lies  in  the  basic  understanding  of  the  bubble 

growth  mechanisms  during  subcooled  nucleate  boiling. 

5.4.1  Validation  of  the  Model  by  the  Comparison  of  Growth 
Rate 

To  validate  the  present  numerical  computation  and  to 
understand  the  effect  of  subcooling  on  the  growth  rate. 
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comparisons  are  made  between  the  predicted  and  the  measured 
R(t)  by  Ellion  (1954)  and  Cochran  and  Aydelott  (1966) . 

Figure  5.4a  shows  comparisons  in  the  high  subcooling  range, 
ATs^jbo”^^*  data  are  taken  from  Ellion  (1954) 

who  used  distilled  degassed  water  as  the  boiling  liquid. 

The  heating  surface  consists  of  a 304  stainless-steel  strip 
with  a thickness  of  0.11  mm.  For  each  subcooling  there  are 
5 sets  of  data  at  the  same  boiling  conditions.  The  data 
shown  in  Figure  5.4a  consist  of  a set  which  falls  within  the 
average  value  of  the  experimental  results.  Close  agreement 
between  the  predicted  and  the  measured  growth  rates  is 

t 

observed  for  most  of  the  time.  Figure  5.4a  also  shows  that 
the  maximum  R(t)  and  the  bubble  lifetime  increase  with  an 
increase  in  bulk  liquid  temperature.  In  the  later  stage  of 
collapse,  the  bubble  radius  is  over-predicted.  It  may  be 
due  to  an  improper  specification  of  the  shape  factor,  C(t) , 
at  later  time.  In  Figure  5.4b  the  computed  bubble  volume  is 
compared  with  the  measured  data  of  Cochran  and  Aydelott 
(1966).  The  subcooling  AT^^j^g  ranges  from  7.4°C~21.4°C,  and 
the  boiling  liquid  is  water.  The  heater  material  is  a 0.13 
mm  thick  chromel  strip.  Good  agreement  for  the  bubble 
volume  is  observed  for  most  of  the  time.  A little 
discrepancy  is  observed  at  the  later  stage  of  collapse.  It 
is  also  evident  that  the  vapor  bubble  has  a larger  volume 
when  the  bulk  liquid  temperature  is  higher.  From  Figures 
5.4a-b  we  see  that  when  the  subcooling  decreases,  the  bubble 
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Figure  5.4a  Comparison  of  the  bubble  growth  rate  between 
the  present  model  and  the  experimental  data  of  Ellion 
(1954)  . 


Volume 
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Figure  5.4b  Comparison  of  the  vapor  bubble  volume  between 
the  present  model  and  the  experimental  data  of  Cochran  and 
Aydelott  (1966). 
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collapse  rate  is  reduced  at  the  later  stage. 

5.4.2  Behavior  of  the  Heat  Transfer  Rate  to  the  Bubble 
during  Growth 

In  Figure  5.5  the  numerical  results  for  the  total  heat 
transfer  rate  between  the  liguid  microlayer  and  the  vapor 
bubble  dome  are  shown.  The  boiling  condition  corresponds  to 
that  of  Ellion  (1954)  with  Tg^jbo=44 . 4 °C.  It  is  evident  that 
the  microlayer  continuously  supplies  energy  to  the  bubble. 

At  the  early  stage,  the  bubble  grows  rapidly  within  the 
thermal  boundary  layer,  in  which  only  a small  portion  of  the 
top  of  the  dome  is  exposed  to  the  subcooled  bulk  liquid. 
Because  of  the  rapid  increase  in  the  bubble  surface  area, 
the  total  heat  transfer  rate  increases  during  the  early 
stage.  In  the  middle  stage  of  the  bubble  growth,  more  of 
the  top  portion  of  the  bubble  surface  area  is  exposed  to  the 
subcooled  liquid  which  results  in  decreasing  the  heat 
transfer  rate  to  the  bubble.  As  the  total  heat  transfer 
rate  becomes  negative,  the  bubble  starts  to  collapse.  When 
the  bubble  is  collapsing  at  the  later  stage,  an  increasing 
portion  of  the  surface  area  is  exposed  to  the  thermal 
boundary  layer  again.  This  results  in  a decreasing  rate  of 
bubble  collapse. 

5.4.3  Thermal  Field  in  the  Solid  Heater  during  the  Growth 
Figure  5.6  shows  the  four  computed  temperature  contours 

at  Ja=83  and  Fo=0.1  for  t=0.1,  0.25,  0.5,  1.  For  Fo=0.1, 
the  heater  is  considered  to  be  thick  and  the  asymptotic 
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Figure  5.5  The  total  heat  transfer  rate  between  the  liquid 
microlayer  and  the  vapor  bubble  dome. 
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Figure  5.6  The  four  computed  temperature  contours  at  Ja=83  and  Fo— 0.1  for 
r=0.1,  0.25,  0.5,  and  1. 
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behavior  (Chapter  3)  is  valid  for  most  of  the  time.  The 
circular  markers  are  the  locations  of  the  microlayer  edge. 
During  the  bubble  growth  period,  it  is  easy  to  see  that  the 
temperature  contours  in  the  solid  penetrate  down  and  move  to 
the  right  in  a similar  way  as  the  bubble  grows.  When  the 
bubble  collapses,  the  temperature  contours  move  to  the  left 
following  the  movement  of  the  microlayer  edge. 

Figure  5.7  shows  the  development  of  the  thermal  field 
at  the  surface  of  the  solid  heater  at  Fo=0.1  and  Ja=83  in 
which  Tg^j^jQ=44 . 4 °C.  The  circular  markers  are  the  locations 
of  the  microlayer  edge.  During  bubble  growth,  the 
temperature  behaves  similarly  to  the  saturated  boiling  case 
(Chapter  3) . It  is  clearly  seen  that  the  thermal  field 
moves  to  the  right  following  the  bubble  growth.  During  the 
bubble  collapse  period,  the  surface  temperature  quickly 
recovers  and  the  thermal  field  shifts  back  to  the  left. 

5.4.4  Thermal  Field  of  the  Bulk  Liquid 

Figure  5.8  shows  the  computed  isotherms  in  the 
subcooled  liquid  surrounding  a growing  and  collapsing  vapor 
bubble.  It  is  seen  that  there  exists  an  unsteady  thermal 
boundary  layer  around  the  vapor  bubble.  The  crease  of  the 
temperature  contours  is  adjacent  to  the  bubble  surface  in  a 
very  thin  thermal  layer.  These  results  are  very  difficult 
to  obtain  experimentally.  Ellion  (1954)  has  stated  that 
there  exists  an  unsteady  thermal  boundary  layer  contiguous 
to  the  vapor  bubble  during  the  bubble  growth  and  collapse. 
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Figure  5.7  Development  of  the  thermal  layer  at  the  surface 

of  the  solid  heater. 
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Figure  5.8  The  computed  isotherms  in  the  subcooled  liquid  surrounding  a growing 
and  collapsing  vapor  bubble  at  four  different  times. 
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Lee  and  Nydahl  (1989)  investigated  the  effect  of 
hydrodynamics  surrounding  the  dome  of  the  growing  bubble  by 
assuming  a constant  wall  temperature  in  saturated  boiling. 
Their  numerical  results  also  showed  the  folding  of  the 
temperature  contours  near  the  bubble  surface.  During  the 
early  stage  of  growth,  we  indeed  see  that  most  of  the  bubble 
surface  is  exposed  to  the  thermal  boundary  layer,  while  a 
portion  of  the  bubble  surface  is  outside  the  thermal 
boundary  layer  when  the  bubble  is  collapsing.  Figure  5.8 
also  showed  that  the  folding  of  the  temperature  contours 
near  the  wall  always  exist  during  bubble  growth  and 
collapse.  The  folding  is  associated  with  the  temperature 
distribution  near  the  bubble  surface  and  will  be  further 
elucidated  by  examining  the  asymptotic  solution. 

5.4.5  Asymptotic  Behavior  of  during  Rapid  Growth 

To  further  understand  the  thermal  boundary  layer  near 
the  bubble  surface,  6^^  based  on  the  asymptotic  solution  and 
numerical  solution  are  compared.  Figure  5.9  shows  the 
comparison  for  dy,  between  the  asymptotic  solutions  and 
numerical  solutions  at  t=0.001,  0.01,  0.1,  and  0.3  for 
\p=0° , 30°,  and  75°.  It  is  clearly  seen  that  the  numerical 
solutions  compare  well  with  the  asymptotic  solutions  at 
T=0.001  and  0.01.  As  mentioned  early,  the  inner  solution  is 
only  valid  when  r-«^l.  At  t=0.1  and  0.3,  R! {t^R)  is  egual  to 
0.3  and  1.2.  Thus,  the  inner  solution  given  by  equation 
(5.36),  is  no  longer  accurate.  However  the  outer  solution 
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Figure  5.9  Comparison  of  the  asymptotic  solutions  and  the 
numerical  solutions  at  t=0.001,  0.01,  0.1,  and  0.3  for 
\p=0° , 30°,  and  75°. 
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should  remain  valid  because  A>1  is  the  only  requirement. 
Figure  5.9  also  shows  the  detailed  structure  of  the  thin 
thermal  boundary  layer  adjacent  to  the  vapor  bubble  surface. 
For  example,  at  t=0.001  and  0.01,  = 14280  and  10950, 

and  the  thermal  boundary  layer  scale  is  5*=A'^^^=0 . 0084  and 
0.0096.  It  is  easily  seen  that  the  thermal  boundary  layer 
thickness  in  Figure  5.9  at  t=0.001  and  0.01  is  the  same 
order  as  5 . It  increases  with  time  because  A decreases 
with  time. 

At  T=0.001  and  0.01,  the  entire  bubble  surface  is 
within  the  thermal  boundary  layer  and  6^^  exhibits  non- 
monotonic behavior  at  i/'=0°  and  ^=30°.  This  non-monotonic 
behavior  is  reflected  by  the  folding  of  the  thermal 
contours.  At  t=0.1  and  0.3,  the  temperature  profile  is 
monotonic  at  ^=0°  and  30°.  The  temperature  profile  is 
always  non-monotonic  at  ^=75°.  This  explains  why  the 
folding  of  the  temperature  contours  near  the  wall  always 
exists  during  bubble  growth  and  collapse. 

5.4.6  Velocity  Field  of  the  Bulk  Liquid 

Figure  5.10  shows  the  computed  velocity  vectors  and  the 
streamlines.  At  t=0.25,  the  flow  field  resembles  that  of  a 
source  flow.  This  source  flow  mainly  results  from  the 
growth  of  the  bubble.  At  t<0.4,  the  flow  field  is  due  to 
the  combination  of  a source  and  a moving  doublet;  the  moving 
doublet  is  caused  by  the  motion  of  the  bubble  center  in  the 
z-direction.  At  t=0.4,  the  velocity  vectors  are  almost 
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Figure  5.10  The  computed  velocity  vectors  and  streamlines 
at  four  different  times. 
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zero.  At  this  point,  the  bubble  growth  rate  is  close  to 
zero,  and  the  bubble  begins  to  collapse.  At  t=0.5,  the  flow 
field  can  be  interpreted  as  a combination  of  a sink  and  a 
doublet  that  is  moving  upwards.  The  doublet  is  moving 
upward  because  dh(t)/dt>0.  At  this  instance  the  vapor 
bubble  continues  to  collapse.  Its  shape  is  becoming  more 
spherical.  Between  r=0.5  and  r=0.75,  the  bubble  collapse 
rate  increases  and  dh(t)/dt  becomes  negative  at  some 
instance.  At  t=0.75,  the  flow  field  is  mainly  due  to  the 
sink  as  the  bubble  continues  to  collapse. 

5 . 5 Summary 

Detailed  theoretical  and  numerical  analyses  are 
presented  to  elucidate  the  effect  of  the  energy  transfer  and 
hydrodynamics  on  the  growth  rate  of  a vapor  bubble  in  highly 
subcooled  pool  boiling.  The  bubble  growth  or  collapse  rates 
are  obtained  by  combining  the  energy  transfered  from  the 
liquid  microlayer  and  from  the  bulk  liquid.  The  microlayer 
wedge  angle,  0,  and  the  relative  growth  rate  (at  t=t„)  , 
are  two  empirically  determined  constants  which  are  used  in 
the  numerical  simulation.  There  only  exists  one  combination 
for  0 and  to  obtain  the  best  fit  for  the  bubble  radius. 

As  long  as  0 and  correctly  specified,  the  predicted 

w 

bubble  growth  rate  compares  well  with  the  reported 
experimental  data.  The  numerical  results  also  reveal 
elegance  of  the  energy  transfer  process  in  the  solid  heater 
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and  the  subcooled  bulk  liquid  during  the  bubble  growth  and 
collapse.  The  rapid  variation  of  the  thermal  field  of  the 
bulk  liquid  during  the  rapid  growth  of  the  bubble  is 
elucidated  by  short  time  asymptotic  solution  that  is 
uniformly  valid  in  the  bulk  liquid.  The  folding  of  the 
temperature  contours  is  explained  via  the  asymptotic  and 
numerical  solution  of  near  the  bubble  surface.  It  is 
shown  that  the  velocity  field  in  the  subcooled  bulk  liquid 
has  a significant  effect  on  the  vapor  bubble  growth  and 
collapse.  The  present  study  for  vapor  bubble  growth  is 
limited  to  pool  boiling  in  the  isolated  bubble  regime. 


CHAPTER  6 


CONCLUSIONS  AND  RECOMMENDATIONS 
FOR  FUTURE  RESEARCH 

This  chapter  summarizes  the  major  research 
accomplishments  from  this  study  and  provides  suggestions  for 
future  work. 


6 . 1 Achievements 

This  research  has  significantly  progressed  the 
understanding  of  the  bubble  growth  process  during  nucleate 
pool  boiling. 

1.  A numerical  analysis  has  been  carried  out  to  study 
vapor  bubble  growth  in  saturated  heterogeneous  pool  boiling. 
The  predicted  bubble  growth  rate  in  saturated  pool  boiling 
compares  very  well  with  the  reported  experimental  data  over 
a wide  range  of  conditions.  The  empirically  determined 
correlation  for  C^  suggests  that  the  liquid  Prandtl  number 
has  little  effect  on  the  bubble  growth  rate.  Thus,  the 
viscosity  of  the  liquid  plays  little  role  in  the  bubble 
growth  process.  A subsequent  estimate  for  the  Reynolds 
number  indicates  that  the  liquid  motion  induced  by  the 
expanding  bubble  is  mainly  of  inviscid  nature. 

2.  The  dimensional  analysis  for  the  governing  equation 
and  boundary  conditions  shows  that  there  exists  four 
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relevant  dimensionless  parameters  for  bubble  growth:  Jacob 
number  Ja,  Fourier  number  Fo,  the  liquid-solid  conductivity 
ratio  K,  and  the  liquid-solid  thermal  diffusivity  ratio  a. 
The  Jacob  number  is  the  most  important  one  affecting  C,  C^, 
and  the  dimensionless  growth  rate  R{z).  The  numerical 
results  in  the  systematic  investigation  assist  in  the  basic 
understanding  of  bubble  growth  mechanisms  in  saturated 
heterogeneous  boiling.  It  is  also  shown  that  the  liquid 
temperature  within  the  microlayer  is  practically  linearly 
during  the  entire  growth  period. 

3.  A simplified  bubble  growth  model  is  developed  in 
saturated  heterogeneous  pool  boiling.  A first  order 
ordinary  differential  equation  is  obtained  for  the  bubble 
growth  rate.  Good  agreement  between  the  predicted  and 
measured  bubble  radius  is  obtained  over  a wide  range  of 
boiling  conditions.  The  present  model  is  a simple  and 
apparently  reliable  alternative  to  the  previous  finite 
difference  computation,  and  provides  modelers  of  boiling 
heat  transfer  with  a very  useful  tool. 

4.  Detailed  theoretical  and  numerical  analyses  have 

been  considered  for  predicting  vapor  bubble  growth  in 
subcooled  heterogeneous  boiling.  Two  empirical  constants, 
the  microlayer  wedge  angle  <p  the  relative  growth  rate, 

fx  , need  to  be  input  in  order  to  start  the  computation.  As 
long  as  the  two  empirical  parameters  are  correctly 
specified,  good  agreement  can  be  obtained  between  the 
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predicted  and  measured  bubble  radius.  It  is  noted  that 
there  is  only  one  combination  for  and  that  provides  the 

w 

best  fit  to  the  bubble  growth  rate  data.  An  asymptotic 
solution  for  small  time  has  been  conducted  for  the  subcooled 
bulk  liquid.  The  numerical  solutions  compare  well  with  the 
asymptotic  solution.  The  computational  study  for  the  solid 
and  the  subcooled  bulk  liquid  thermal  fields  elucidates  the 
energy  transfer  process  during  the  bubble  growth  and 
collapse.  It  is  revealed  that  the  subcooled  bulk  liquid 
flow  has  a significant  effect  on  the  vapor  bubble  growth  and 
collapse  process. 


6 . 2 Suggestions  for  Future  Work 

The  motivation  for  the  present  research  is  to  gain  a 
fundamental  understanding  of  the  vapor  bubble  growth  process 
in  saturated  and  subcooled  nucleate  pool  boiling.  Although 
significant  progress  has  been  made  toward  achieving  goal, 
uncertainties  remain.  In  order  to  build  on  the  present 
research  and  continue  further  studies,  the  following  future 
investigations  are  recommended: 

1.  Develop  a reliable  model  for  predicting  vapor 
bubble  incipience. 

2.  Extend  the  finite  difference  computation  to  include 
the  simultaneous  bubble  growth  of  multiple  bubbles  in 
heterogeneous  pool  boiling. 

3.  Develop  a reliable  bubble  departure  model  in 
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subcooled  heterogeneous  pool  boiling. 

4.  Extend  the  bubble  growth  model  to  include  flow 
boiling.  In  the  case  of  flow  boiling,  the  liquid  inertia 
will  distort  the  bubble  shape.  In  addition  the  bubble 
slides  along  the  heating  surface  prior  to  detachment. 
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